Modifying fs, fmin, and fmax for DirectedTransferFunction

[withercp]

Hi! I have a question about the DirectedTransferFunction implementation. I am using intracranial EEG data, and I would like to compute DTF on a 5 second epoch sampled at 500Hz. My data has shape (60,2500), where 60 is my # of channels and 500*5 is the # of timepoints. I would like to compute DTF in the beta band specifically (13-30Hz). Here is how I have formatted my config:

.statistics.spectral:
  DirectedTransferFunction:
    labels:
      - directed
      - linear
      - unsigned
      - bivariate
      - frequency-dependent
    dependencies:
    configs:
      - fs: 500
        fmin: 13
        fmax: 30
        statistic: mean

Does this work properly or has it not been implemented yet? Additionally, is it possible to get out the full spectrum DTF or do I have to compute the mean/max?

Thank you so much!

[joshuabmoore]

Hi @withercp,

If your goal is to compute the directed transfer function in the beta band, then setting the sampling frequency, min, and max frequencies like this is correct for your specific problem. You might get the warning Multiple sampling frequencies not yet handled, however this should not prevent the computation of the SPI.

While Pyspi is designed to only return a single scalar value (summary statistic), the results from the computation are stored internally and therefore you can access the full DTF spectrum from the calculator as follows:

Run the calculator compute method as usual with your configuration specified above:

calc = Calculator(dataset=<your_data>, configfile=<your_config>.yaml)
calc.compute()

Now, you can access the frequencies and DTF with:

freq = calc.dataset.spectral_bv['freq']
dtf_i_j = calc.dataset.spectral_bv[("directed_transfer_function", <i>, <j>)]

Replace i and j with the processes of interest (which in your case would be the EEG channels). For example, to extract the DTF computed from source channel 0 to target channel 1, use:

calc.dataset.spectral_bv[("directed_transfer_function", 0, 1)]

The result will have shape (n_windows, n_fft_samples, n_signals, n_signals) as per the documentation here.

You can then go ahead and plot the spectrum:

import matplotlib.pyplot as plt
fig, ax = plt.subplots()
ax.plot(freq, dtf_i_j [0, :, 0, 1]) # source 0 to target 1
ax.set_xlabel('freq [Hz]')
ax.set_ylabel('DTF')
plt.show()