Suggestions on preparing input for physio_calc.py

Hello AFNI experts,

I would like to use physio_calc.py to prepare a 1D file for use with -ricor_regs. My understanding of the documentation is that specifying -dset_epi helps minimize human error by extracting relevant information directly from the fMRI data. However, my physiological measurement files, acquired on a Siemens Terra 7T scanner, contain data from all runs of an entire single acquisition, including respiration (*.resp), cardiac pulsation (*.puls), and the trigger log file (*.ext). I am wondering whether physio_calc.py can automatically segment the data into individual runs. In other words, giving a list of fMRI runs to -dset_epi as in afni_proc.py. If not, do you know of any other tools that can perform this task? I am thinking to generate *_physio.tsv files for each run following the BIDS convention and then create a shell script to loop through each run, applying physio_calc.py to the files (*_physio.tsv and *.nii) in the func folder. Does this approach make sense to you? Thank you very much.

One potential "mistake" in the physio_calc.py manual is that:

It appears that the number is inserted into the series, in which case,
5000 values could simply be removed rather than replaced by an
interpolation of the two adjacent values, using the option
'remove_val_list ..'.

This is also what I learned, which means 5000 should be removed from the data

-extra_fix_list FVAL1 [FVAL2 ...]
List of one or more values that will also be
considered 'bad' if they appear in the physio time
series, and replaced with interpolated values

but in the example at the bottom, '-extra_fix_list 5000' is used, which will interpolate the trace that may create a considerable amount of time drifting (extra data points) since there are actually many 5000s. I guess 'remove_val_list' should be used, or usually just assume 5000s are already removed when praparing the input for physio_calc.py. I guess this will not cause issues in most of the cases, as people who prepare the physio files usually remove 5000.

After searching for a while, I couldn’t find a tool specifically designed for Siemens physiological recordings (I don't want to install dependences for this task). However, this thread provides enough information to create a simple script that automatically extracts physiological data chunks from a file containing the resp and pulse of the entire acquisition session. The key is to compute the timestamp for each data point using LogStartMDHTime and LogStopMDHTime found in the footer of the *.puls and *.resp files separately (as they may not always have the exact same sampling rate). Then, the appropriate chunk can be selected based on the AcquisitionTime tag in each BOLD JSON file, TR and number of volumes.

If a *.ext file (containing TR triggers) is available, the process becomes much easier since the data points in the three physiological files should be "well-aligned" when working backward from the last recorded data point. This approach assumes a sufficiently high sampling frequency, which, in our case, is 400 Hz for *.ext and 399.9x Hz for *.puls and *.resp, fast enough to ignore minor time differences between aligned data points.

It seems for me aligning the data points backward from the last recorded sample appears to be more reliable. This is because the LogStopMDHTime values in all three files are identical, whereas the LogStartMDHTime values differ a bit in my case.

Yes, I believe the 5000 should be removed from the input. Based on the documentation here, there are several special values - 5000, 5002, 5003, 6000.

https://github.com/translationalneuromodeling/tapas/blob/master/PhysIO/code/readin/tapas_physio_siemens_line2table.m

1 Like

Hi Daniel (@dglen), I get an error ** ERROR: inconsistent slice times entered, it may be related to the slice timing option I used. Could you please help me to understand this error?

I did the following

physio_calc.py                                                           \
    -phys_file           '${path-physiofile}_run-1_physio.tsv'                               \
    -phys_json               '${path-physio}_physio.json'                                             \
    -dset_epi            '${path-dset}_run-1_echo-1_bold.nii.gz'                                        \
    -dset_slice_pattern  alt+z                                           \
    -do_fix_nan                                                          \
    -out_dir             '/mnt/d/'                                         \
    -prefix              'physio'

and got an error

++ Slice pattern from cmd line: 'alt+z'
+* Difference in list elements:
A[1] = 1.0133333333333334
B[1] = 1.04
+* Difference in list elements:
A[2] = 0.02666666666666667
B[2] = 0.08
+* Difference in list elements:
A[3] = 1.04
B[3] = 1.12
+* Difference in list elements:
A[4] = 0.05333333333333334
B[4] = 0.16
+* Difference in list elements:
A[5] = 1.0666666666666667
B[5] = 1.2
+* Difference in list elements:
A[6] = 0.08
B[6] = 0.24
+* Difference in list elements:
A[7] = 1.0933333333333333
B[7] = 1.28
+* Difference in list elements:
A[8] = 0.10666666666666667
B[8] = 0.32
+* Difference in list elements:
A[9] = 1.12
B[9] = 1.36
+* Difference in list elements:
A[10] = 0.13333333333333333
B[10] = 0.4
+* Difference in list elements:
A[11] = 1.1466666666666667
B[11] = 1.4375
+* Difference in list elements:
A[12] = 0.16
B[12] = 0.48
+* Difference in list elements:
A[13] = 1.1733333333333333
B[13] = 1.5175
+* Difference in list elements:
A[14] = 0.18666666666666668
B[14] = 0.56
+* Difference in list elements:
A[15] = 1.2
B[15] = 1.5975
+* Difference in list elements:
A[16] = 0.21333333333333335
B[16] = 0.64
+* Difference in list elements:
A[17] = 1.2266666666666666
B[17] = 1.6775
+* Difference in list elements:
A[18] = 0.24
B[18] = 0.72
+* Difference in list elements:
A[19] = 1.2533333333333334
B[19] = 1.7575
+* Difference in list elements:
A[20] = 0.26666666666666666
B[20] = 0.8
+* Difference in list elements:
A[21] = 1.28
B[21] = 1.8375
+* Difference in list elements:
A[22] = 0.29333333333333333
B[22] = 0.88
+* Difference in list elements:
A[23] = 1.3066666666666666
B[23] = 1.9175
+* Difference in list elements:
A[24] = 0.32
B[24] = 0.96
+* Difference in list elements:
A[25] = 1.3333333333333333
B[25] = 0.0
+* Difference in list elements:
A[26] = 0.3466666666666667
B[26] = 1.04
+* Difference in list elements:
A[27] = 1.36
B[27] = 0.08
+* Difference in list elements:
A[28] = 0.37333333333333335
B[28] = 1.12
+* Difference in list elements:
A[29] = 1.3866666666666667
B[29] = 0.16
+* Difference in list elements:
A[30] = 0.4
B[30] = 1.2
+* Difference in list elements:
A[31] = 1.4133333333333333
B[31] = 0.24
+* Difference in list elements:
A[32] = 0.4266666666666667
B[32] = 1.28
+* Difference in list elements:
A[33] = 1.44
B[33] = 0.32
+* Difference in list elements:
A[34] = 0.4533333333333333
B[34] = 1.36
+* Difference in list elements:
A[35] = 1.4666666666666666
B[35] = 0.4
+* Difference in list elements:
A[36] = 0.48
B[36] = 1.4375
+* Difference in list elements:
A[37] = 1.4933333333333334
B[37] = 0.48
+* Difference in list elements:
A[38] = 0.5066666666666667
B[38] = 1.5175
+* Difference in list elements:
A[39] = 1.52
B[39] = 0.56
+* Difference in list elements:
A[40] = 0.5333333333333333
B[40] = 1.5975
+* Difference in list elements:
A[41] = 1.5466666666666666
B[41] = 0.64
+* Difference in list elements:
A[42] = 0.56
B[42] = 1.6775
+* Difference in list elements:
A[43] = 1.5733333333333333
B[43] = 0.72
+* Difference in list elements:
A[44] = 0.5866666666666667
B[44] = 1.7575
+* Difference in list elements:
A[45] = 1.6
B[45] = 0.8
+* Difference in list elements:
A[46] = 0.6133333333333333
B[46] = 1.8375
+* Difference in list elements:
A[47] = 1.6266666666666667
B[47] = 0.88
+* Difference in list elements:
A[48] = 0.64
B[48] = 1.9175
+* Difference in list elements:
A[49] = 1.6533333333333333
B[49] = 0.96
+* Difference in list elements:
A[50] = 0.6666666666666666
B[50] = 0.0
+* Difference in list elements:
A[51] = 1.68
B[51] = 1.04
+* Difference in list elements:
A[52] = 0.6933333333333334
B[52] = 0.08
+* Difference in list elements:
A[53] = 1.7066666666666668
B[53] = 1.12
+* Difference in list elements:
A[54] = 0.72
B[54] = 0.16
+* Difference in list elements:
A[55] = 1.7333333333333334
B[55] = 1.2
+* Difference in list elements:
A[56] = 0.7466666666666667
B[56] = 0.24
+* Difference in list elements:
A[57] = 1.76
B[57] = 1.28
+* Difference in list elements:
A[58] = 0.7733333333333333
B[58] = 0.32
+* Difference in list elements:
A[59] = 1.7866666666666666
B[59] = 1.36
+* Difference in list elements:
A[60] = 0.8
B[60] = 0.4
+* Difference in list elements:
A[61] = 1.8133333333333332
B[61] = 1.4375
+* Difference in list elements:
A[62] = 0.8266666666666667
B[62] = 0.48
+* Difference in list elements:
A[63] = 1.84
B[63] = 1.5175
+* Difference in list elements:
A[64] = 0.8533333333333334
B[64] = 0.56
+* Difference in list elements:
A[65] = 1.8666666666666667
B[65] = 1.5975
+* Difference in list elements:
A[66] = 0.88
B[66] = 0.64
+* Difference in list elements:
A[67] = 1.8933333333333333
B[67] = 1.6775
+* Difference in list elements:
A[68] = 0.9066666666666666
B[68] = 0.72
+* Difference in list elements:
A[69] = 1.92
B[69] = 1.7575
+* Difference in list elements:
A[70] = 0.9333333333333333
B[70] = 0.8
+* Difference in list elements:
A[71] = 1.9466666666666668
B[71] = 1.8375
+* Difference in list elements:
A[72] = 0.96
B[72] = 0.88
+* Difference in list elements:
A[73] = 1.9733333333333334
B[73] = 1.9175
+* Difference in list elements:
A[74] = 0.9866666666666667
B[74] = 0.96
** ERROR: inconsistent slice times entered

However, if I remove -dset_slice_pattern, it runs without any error, as the dset contain slice_timing added via abids_tool.py -add_slice_times -input "$sub_dir"/func/*bold.nii.gz right after DICOM to BIDS converstion. The slice timing is the following, and I think it is alt+z.

"dset_slice_times": [
        0.0,
        1.04,
        0.08,
        1.12,
        0.16,
        1.2,
        0.24,
        1.28,
        0.32,
        1.36,
        0.4,
        1.4375,
        0.48,
        1.5175,
        0.56,
        1.5975,
        0.64,
        1.6775,
        0.72,
        1.7575,
        0.8,
        1.8375,
        0.88,
        1.9175,
        0.96,
        0.0,
        1.04,
        0.08,
        1.12,
        0.16,
        1.2,
        0.24,
        1.28,
        0.32,
        1.36,
        0.4,
        1.4375,
        0.48,
        1.5175,
        0.56,
        1.5975,
        0.64,
        1.6775,
        0.72,
        1.7575,
        0.8,
        1.8375,
        0.88,
        1.9175,
        0.96,
        0.0,
        1.04,
        0.08,
        1.12,
        0.16,
        1.2,
        0.24,
        1.28,
        0.32,
        1.36,
        0.4,
        1.4375,
        0.48,
        1.5175,
        0.56,
        1.5975,
        0.64,
        1.6775,
        0.72,
        1.7575,
        0.8,
        1.8375,
        0.88,
        1.9175,
        0.96
    ]

Although physio_calc.py runs without specifying the slice pattern, I would prefer to understand this error. Physiological data are sort of "stable" oscillations, making it difficult to detect potential errors, especially when the issue involves wrong timing (hard to tell an error from the trace plot).

Another question: I have multi-echo data, and the slice timing for each echo .nii image is the same. Should I include the echo time in the computation to accurately determine the timing for each slice of every echo? Is there any example from AFNI? I am thinking to create .json file for each echo image, and put the echo time for the StartTime, not sure if this is correct. In the afni_proc.py examples, Example 13. Complicated ME, surface-based resting state example. used multi-echo data with -ricor_regs_nfirst 2 with ricor, does this automatically account for echo timing? if so, any trick to prepare 1D file for -ricor_regs? should the timing be aligned to the 1st echo?

{
  "StartTime": 0,
  "SamplingFrequency": 400,
  "Columns": [
    "cardiac",
    "respiratory"
  ]
}

Thank you very much

Those slice times repeat. Is this data also multi-band=2?

1 Like

Thanks, these are multiband 3, here are some relevent JSON tags. I just did DICOM to BIDS and then abids_tool.py -add_slice_times -input "$sub_dir"/func/*bold.nii.gz which add slice timing from the JSON file to nifti history.

   "SAR":0.128952,
    "EchoNumber":2,
    "EchoTime":0.02965,
    "RepetitionTime":2,
    "FlipAngle":67,
    "PartialFourier":0.75,
    "BaseResolution":118,
    "ShimSetting":[
        151,
        -83,
        23,
        -78,
        -9,
        -130,
        -127,
        -271
    ],
    "TxRefAmp":240,
    "PhaseResolution":1,
    "ReceiveCoilName":"1Tx32Rx_Head",
    "CoilString":"A32",
    "PulseSequenceDetails":"%CustomerSeq%\\cmrr_mbep2d_bold",
    "WipMemBlock":"e5eba07d-9c0c-4db3-b185-6444cd4a05ea||Sequence: R016 ve12u/SP01_R016a r/af3652dd; Apr  1 2021 14:05:17 by eja",
    "RefLinesPE":54,
    "CoilCombinationMethod":"Sum of Squares",
    "ConsistencyInfo":"N4_VE12U_LATEST_20181126",
    "MatrixCoilMode":"GRAPPA",
    "MultibandAccelerationFactor":3,
    "PercentPhaseFOV":100,
    "PercentSampling":100,
    "EchoTrainLength":29,
    "PhaseEncodingSteps":88,
    "AcquisitionMatrixPE":118,
    "ReconMatrixPE":118,
    "BandwidthPerPixelPhaseEncode":43.834,
    "ParallelReductionFactorInPlane":3,
    "EffectiveEchoSpacing":0.000193333,
    "DerivedVendorReportedEchoSpacing":0.00058,
    "TotalReadoutTime":0.02262,
    "PixelBandwidth":2355,
    "DwellTime":1.8e-06,
    "PhaseEncodingDirection":"j-",
    "SliceTiming":[
        0,
        1.04,
        0.08,
        1.12,
        0.16,
        1.2,
        0.24,
        1.28,
        0.32,
        1.36,
        0.4,
        1.44,
        0.48,
        1.52,
        0.56,
        1.6,
        0.64,
        1.68,
        0.72,
        1.76,
        0.8,
        1.84,
        0.88,
        1.92,
        0.96,
        0,
        1.04,
        0.08,
        1.12,
        0.16,
        1.2,
        0.24,
        1.28,
        0.32,
        1.36,
        0.4,
        1.44,
        0.48,
        1.52,
        0.56,
        1.6,
        0.64,
        1.68,
        0.72,
        1.76,
        0.8,
        1.84,
        0.88,
        1.92,
        0.96,
        0,
        1.04,
        0.08,
        1.12,
        0.16,
        1.2,
        0.24,
        1.28,
        0.32,
        1.36,
        0.4,
        1.44,
        0.48,
        1.52,
        0.56,
        1.6,
        0.64,
        1.68,
        0.72,
        1.76,
        0.8,
        1.84,
        0.88,
        1.92,
        0.96
    ],
    "ImageOrientationPatientDICOM":[
        1,
        0,
        0,
        0,
        0.99863,
        -0.052336
    ],
    "ImageOrientationText":"Tra>Cor(-3.0)",
    "InPlanePhaseEncodingDirectionDICOM":"COL",

You can't use traditional slice timing pattern names like "alt+Z" with multi-band data. The patterns can be stored in the EPI files as a pattern with "3drefit -Tslices", applied in "3dTshift -tpattern" or gotten from the json file. NIFTI doesn't support the slice times directly in the NIFTI header, only as an extension. I don't think we read that extension yet, but physio_calc should read the slice times from the json file, I believe.

1 Like

Yes, abids_tool.py -add_slice_times adds the slice timing to nifti history which will be read for slice timing correction and physio_calc.

In the example 13 copied below, should I use "StartTime": 0, for physio.json used in physio_calc? assuming afni_proc.py will consider the echo time when doing ricor. Or should I using the echo time of the 1st echo for "StartTime"? It seems the .slibase.FT.r?.1D file is one column per slice per regressor, but which echo this file is referring to?

afni_proc.py                                                         \
    -subj_id                   FT.complicated                        \
    -dsets_me_echo             FT/FT_epi_r?+orig.HEAD                \
    -dsets_me_echo             FT/FT_epi_r?+orig.HEAD                \
    -dsets_me_echo             FT/FT_epi_r?+orig.HEAD                \
    -echo_times                11 22.72 34.44                        \
    -blip_forward_dset         'FT/FT_epi_r1+orig.HEAD[0]'           \
    -blip_reverse_dset         'FT/FT_epi_r1+orig.HEAD[0]'           \
    -copy_anat                 FT/FT_anat+orig                       \
    -anat_follower_ROI         FSvent epi FT/SUMA/FT_vent.nii        \
    -anat_follower_erode       FSvent                                \
    -blocks                    despike ricor tshift align volreg     \
                               mask combine surf blur scale regress  \
    -radial_correlate_blocks   tcat volreg                           \
    -tcat_remove_first_trs     2                                     \
    -ricor_regs_nfirst         2                                     \
    -ricor_regs                FT/fake.slibase.FT.r?.1D              \
    -ricor_regress_method      per-run                               \

Interesting question. We don't use the echo times at all in the slice timing calculation. This would typically be a very small fraction of the TR, and every slice will have range of time determined by the echo times. In our example multi-echo demo data (installed with @Install_APMULTI_Demo1_rest ) , the slice time offsets start at 0 and goes according to an alt+Z type pattern with no multi-band.

1 Like

Thanks @dglen, as my slice time interval is 80 ms, which is considerably larger than the TE difference across echoes, and considering that even the fastest physio measurement - pulse trace is not as fast as a few tens of ms, I guess I can assume that the same physio traces can be safely applied to multiecho data.