Tensorpac: Pac values are nan for low frequencies plus weird results

đź’¬ Support & Discussions
1

Hi,

I have applied a phase amplitude coupling analysis to my data. I chose the range 1 to 8 hz for the phase (width=2), and 8 to 120 for the amplitude frequencies (width=4).

When I run the analysis, I am faced with two problems.

  1. for a couple of subjects I get missing data in the comodulogram:

image

I cannot understand why I get these empty columns for low frequencies only for a few subjects, when I applied to all participants the same exact pre-processing steps.

  1. For those participants with full comodulogram (e.g., see image below), I also feel the pattern of results is a bit at odds compared to what I have observed in other studies. Apart from the 0 PAC values at 50Hz (which is noise due to a nearby electronic device - I still have to use a notch filter), it seems strange to me to see these relatively high PAC at lower phase frequencies almost along the whole column of amplitude frequencies. Is there something wrong? For those expert in PAC analysis, can you spot a potential problem?

image

Here is part of the code:

import tensorpac
from tensorpac import Pac
time_delay = (1.5,3)
time_delay_idx = slice(*[np.argmin(abs(phase['1'].times - t)) + i for t,i in zip(time_delay, (0,1))])
sf = 256
pac_ato = []
pac_atb = []
subjects = list(range(1,37))
for i in subjects:
    print(i)
    p_obj = Pac(idpac=(2, 0, 0), f_pha=(1, 8,2,0.2), f_amp=(8,120,4,1))
    pac_electrodes_ato = []
    pac_electrodes_atb = []
    for e in range(0,64): #electrodes
        # condition 1
        pha_p_ato = p_obj.filter(sf, EEG[f'{i}'][act].copy()._data[:,e,:], ftype='phase')
        amp_p_ato = p_obj.filter(sf, EEG[f'{i}'][act].copy()._data[:,e,:], ftype='amplitude')
        pha_delay_ato, amp_delay_ato = pha_p_ato[..., time_delay_idx], amp_p_ato[..., time_delay_idx]
        pac_delay_ato = p_obj.fit(pha_delay_ato, amp_delay_ato).mean(-1)
        pac_electrodes_ato.append(pac_delay_ato)
        # condition 2
        pha_p_atb = p_obj.filter(sf, EEG[f'{i}'][cont].copy()._data[:,e,:], ftype='phase')
        amp_p_atb = p_obj.filter(sf, EEG[f'{i}'][cont].copy()._data[:,e,:], ftype='amplitude')
        pha_delay_atb, amp_delay_atb = pha_p_atb[..., time_delay_idx], amp_p_atb[..., time_delay_idx]
        pac_delay_atb = p_obj.fit(pha_delay_atb, amp_delay_atb).mean(-1)
        pac_electrodes_atb.append(pac_delay_atb)
        
    pac_ato.append(pac_electrodes_ato)
    pac_atb.append(pac_electrodes_atb)

Then to show the graph for one participant I used e.g.:

p_obj.comodulogram(np.mean(np.array(pac_ato)[0],axis=0))

Many thanks for your help!

2

Hi @Catsine

Your comodulogram does seem odd. The coupling between the lowest driver frequency and the carrier frequencies looks like it could be an edge effect.

Sorry to be unhelpful but Tensorpac is not a tool within the MNE ecosystem so I’m not sure anyone here will have any insights into your problem. But as I recall tensorpac does have its own user forum.

You could also try pactools to see if your issues replicate with different software?

3

Hi @scott-huberty,

thanks for your reply. I obtain the same results also with concatenated data, for which there should not be edge effects (I also drop two segments of the data at the beginning and end of the time window of interest). I will try to use their forum, and take a look at the pactools package.

Many thanks for your comment!

1 Like
4

Sorry I couldn’t be of more help! Maybe someone else can. I’m not sure how active tensorpac is these days, trying to replicate the problem with pactools is probably where I’d start!

5

If I get the same weird results with pactools, would you be able to help me? As it seems that I get similar (weird) results… :sweat_smile:

1 Like
6

It’s been a long time since I’ve used either pactools or tensorpac but if you get similar results with pactools and share your code and a file (or if your file is confidential, then if you can replicate the code with one of MNE’s example files), then I can try to take a look!

1 Like
7

Thanks! I will upload something asap

1 Like
8

Hi @scott-huberty,

here I am. I will share with you my Jupyter Notebook code and the dataset (both behaviour and EEG). I included all the information tin the script file. Could you send me your email? I have created a compressed folder. Thanks!

9

Or you can use this link https://drive.google.com/drive/folders/1tGlbgWuoczGI4nTQMyqVjVX0-ruBRIyo?usp=drive_link

10

Hi @scott-huberty,

Did you manage to see the folder I sent you? Is there any problem with the data?

Thanks!

11

Hi @Catsine sorry we had a long weekend in North America so I’ve just taken a look at this.

The code is quite imbricated, which is not meant to be an insult, I know sometimes this is unavoidable for some datasets, but it makes it a lot easier to help if you share a file and a notebook that “works out of the box”, i.e. someone can just run it without having to adjust code, filepaths etc.

You import a lot of libraries that you don’t seem to use (pygazeanalyzer?) and include a hardcoded path to a code base (surface_laplacian) that I don’t have access to. Unfortunately I won’t have time to dig through it to get to the root of the problem.

But below are a couple of naive comments on the code, if they might be helpful:

1.Double check the way you are indexing your epochs:

You are indexing your epochs using the indices of a pandas DataFrame (i.e. epochs[df.index]). I think the assumption is that each row of the dataframe represents behavioural data for 1 epoch in your EEG data. But the number of rows in the dataframe and the number of epochs do not match. So maybe double check that this indexing is correct.

2. Double check the length of the data that you are calculating PAC with

It looks like you are cropping your epochs, and then you further slice the times dimension, before calculating PAC. It makes me wonder if you are calculating PAC on very short chunks of data. Typically you want to have long epochs (i.e. 4-seconds or longer), since Delta (which is your driver frequency) only oscillates 1-4 times per second (depending on how you define its frequency range).

3. Try to simplify your code.

I’d try to simplify the code as much as possible, before calculating PAC. for example, you could drop the surface_laplacian functions from your routine, or extract your epochs without further slicing them, and see if that improves the output of your comodulogram. Then, you can gradually start to re-integrate the code you took out.

Alternatively, caclulating PAC with pactools could help you figure out if your comodulogram issues are a bug in tensorpac or a problem with your preprocessing.

If you are still having trouble, maybe try popping by the next MNE office hours!

1 Like