Dear all,

after using 3dPeriodogram for a longer time now, I stumbled across a strange “problem”. I am interested in computing power spectra and additionally I am interested in computing the slope of a linear regression in the log-frequency log-power transformation of the frequency-domain, namely the power-law exponent (PLE).

Here are two options how to accomplish this, while option 1 is clearly prefered.

Option 1:
- Use 3dPeriodogram on every voxel in AFNI on the preprocessed time-series.
- Cut the 3dPeriodogram results to the desired frequency range via 3dTcat.
- (Optional: use 3dcalc to log the data, or log them later via numpy in Python).
- Extract the results via 3dmaskave into a .1D or .txt file (this provides the average of all voxels’ Periodograms for a desired ROI).
- Load this 3dmaskave result into Python and compute the slope.

Option 2:
- Extract the preprocessed time-series with 3dmaskave choosing the desired ROI.
- Load the .1D or .txt file (the average across all voxels’ time-series) into Python.
- Use scipy Periodogram to transform the time-series into the frequency-domain.
- Log-log the Periodogram (the frequency and data) and compute the slope of log-frequency and log-power.

Both results should yield more or less the same power-spectra and PLE (I heard so). A lab colleague who is very experienced with EEG recommended option 1, since they average the power spectra across the single electrodes as well.
Option 1 yields a PLE of ~ -0.8, while option 2 yields a PLE of roughly twice the value, namely ~ -1.6. This “problem” is not related to my Python script (I am 99% certain about this!), not to a specific subject. But the problem must stem from something that 3dPeriodogram simply does different than scipy Periodogram and scipy Welch.

True, in option 1, we average thousands over thousands of Periodograms (across the voxels, it is a huge ROI), while in option 2 we average thousands over thousands of time-series. But this should not result into such a massive power spectra (see attached images) and PLE difference.

Furthermore, when I take the 3dPeriodogram or 1dFFT results ^2 (square them) in Python before applying log-log and the linear regression on them, they yield approximately the same results as scipy Periodogram or scipy Welch. Is precisely THIS related to the “problem” I am having here? I no longer know which power spectra and PLE values to trust, either the one from option 1 or option 2, since I have no ground truth.

I know that this is probably not directly an AFNI related problem, but the main question here is if 3dPeriodogram and 1dFFT do something so differently compared to the scipy Periodogram method that the power distribution between both methods looks so different (see the attached images).
The options in scipy Periodogram are all matched to 3dPeriodgram (window = Hamming; scaling = Density, nfft = the same, fs = matches to TR, detrend = linear, etc). But I tried out several options anyway, the basic problem remains.

Philipp



Edited 1 time(s). Last edit at 05/27/2022 10:39AM by Philipp.