Matplotlib seems to be listed as an optional requirement but any import of tensorpac immediately imports matplotlib as well.

import tensorpac leads to:

Can it be changed so that importing any module doesn't cause utils.py to run or if matplotlib was a late import so it's only needed when using the plot functions?

Hello, I would like to know whether the resting state EEG needs to be surrogate, and the range of PAC calculated by several methods in this example is not consistent, that is, the color bar on the right side of each subplot, and the analysis of the two groups The difference of PAC on a certain electrode of a person, how to do a cluster-based permutation test, thank you very much

There are no github tags for recent releases. These are very helpful for people who want the complete source rather than just what is included in the pypi sdists.

Hi Etienne, in the README example I tested the Modulation Index method (among others) and all is good, until I remove the noise. Setting the noise to 0 and by having a completely clean and highly-coupled signal breaks the plot results. Is this an expected artifact of calculating values from signals that are simply not there (filtering outside the 10Hz/100Hz does not yield anything), or is there something else that I am overlooking?

Thanks in advance, great job with the toolbox!

Hello,

I just to want to express my appreciation for the amazing work you have done in this project. Truly the level of attention to details is astounding.
I have a question about the code for spectral analysis using Morlet wavelets.
I am seeing computing of the sine wave, as well the Gaussian envelope. I am guessing 'a' is the scaling factor, together with the normalization factor which is added for the wavelet to have a unit energy. I think you are using (1/sqrt(sigma)) as a normalization term which is a popular normalization scheme for morlet wavelets.

Having used this function, although the power spectrum matches what I got from other libraries graphically ( pycwt or matlab old-version toolbox), but the raw values don't seem to be normalized to the [0, 1] range. I am constantly getting values in the range of tens. Comparing both libraries, values align by a factor of 10 to the power of 3, which is oddly enough is the sampling frequency of the signals.

I am aware it is not relevant if you are doing a baseline normalization, but for multiple reasons, I am not able to do that in my current analysis. Can you point to the right direction?

Hi Etienne, in pac.idpac you credit true pac method 5 (Phase-Locking-Value; tensorpac.methods.phase_locking_value()) to Both Lachaux 1999 and W.D. Penny 2008, and in your paper your figure for method 5 credits only W.D. Penny 2008. In their 2008 paper, Penny et al refer to their method as a General Linear Model (GLM) method, which is derived from the envelop-to-signal-correlation (ESC) method. And in fact they compare their GLM method against a PLV method. If your implementation of p.idpac method 5 is based Penny's paper, would it be fair to refer to it as GLM PAC? Alternatively, is the method actually only based on Lacheaux's 1999 paper, which I believe does use a Phase-Locking-Value?

Best,

Scott

hi @EtienneCmb

after talking with Karim this week I took 15 mins to look at your package.

a quick feedback:

• use sphinx-gallery to build the doc with examples
• make sure legend + titles are readable
• your docstrings don't follow numpydoc conventions
• you use sine functions for simulations. It's too simple
• don't use CamelCase but rather snake_case for function names in python
• your code coverage is weak (only 24%)
• add a flake8 test for your codebase. There are a lot of style issues.

you should also talk to @TomDLT working with me as he has started a similar
effort at : https://pactools.github.io/

I have a bunch of animals I sent through tensorpac (hilbert instead of wavelet), and all of their data are nearly the same length, but strangely enough, some took 10 minutes for tensorpac to crunch, some took hours or days. It turns out to boil down to an annoying feature of the hiblert method. If time is not a Hamming number, it takes 3000x longer to run. So here's my fix in spectral.py:

def hilbert_fast(x, axis): 
    ''' If x not factorizable by (2,3,5), takes like 3000x longer '''
     from scipy import fftpack
     if isinstance(x, (list,tuple)):
         x = np.array(x) # concat iterables
     xd = hilbert(x, axis=axis, N=fftpack.next_fast_len(x.shape[axis])) # pad to a hamming number
     xd = xd.swapaxes(0,axis) # swap such that time along 0-axis
     xd = xd[:x.shape[axis]].swapaxes(0,axis) # cut out the non-padded portion and return the axes to their positions
    return  xd

And then you sub the hilbert call with the following:
xd = hilbert_fast(xf, axis=axis + 1) if stype is not None else np.array(xf)

Night and day speed difference for several of my animals.

Whether the color bars are t-values？

Hello,

I was interested in comparing how low vs high PAC coupling looks, and so was interested in the included function to generate signals. However, I noticed that even when assigning high values of coupling (chi =0.0), pvalues aren't detected as significant above the permutations, though this isn't true when generating the signal with wavelet.

For instance, using:

f_pha = 6 # frequency phase for the coupling
f_amp = 90 # frequency amplitude for the coupling
n_epochs = 1 # number of trials

sf = 1000. # sampling frequency
n_times = 10*sf # number of time points

data, time = pac_signals_wavelet(f_pha=f_pha, f_amp=f_amp, noise=0.8,
n_epochs=n_epochs, n_times=n_times, sf=sf, pp=0.0)
p = Pac(idpac=(2, 3, 0), f_pha=(2, 29, 1, .2), f_amp=(30, 120, 5, 1))
xpac = p.filterfit(sf, data2, n_perm=n_perm, n_jobs=4, random_state=0 ).squeeze()
xpac_smean = xpac[pval < .05].mean()
p.comodulogram(xpac, title='Modulation Index (MI)', cmap='viridis', vmin=0., pvalues=pval,
levels=.05)

I am able to produce the figure as given in https://etiennecmb.github.io/tensorpac/auto_examples/stats/plot_pvalues.html#sphx-glr-auto-examples-stats-plot-pvalues-py

However, when generating the signal using tort:
data, time = pac_signals_tort(f_pha=f_pha, f_amp=f_amp, sf=sf, n_times=n_times, n_epochs=n_epochs,
chi=0.0, noise=2.4, dpha=0.0, damp=0.0, rnd_state=0)
p = Pac(idpac=(2, 3, 0), f_pha=(2, 29, 1, .2), f_amp=(30, 120, 5, 1))
xpac = p.filterfit(sf, data, n_perm=n_perm, n_jobs=4,
random_state=0).squeeze()

pval = p.infer_pvalues(p=0.05)
xpac_smean = xpac[pval < .05].mean()
p.comodulogram(xpac, title='Modulation Index (MI) - wavelet', cmap='viridis', vmin=0., pvalues=pval,
levels=.05)

I get a similar region with high MI but it is no longer significant

And this is true with different noise levels, using dcomplex='hilbert' or 'wavelet', different methods of permutation calculation, different numbers of permutations, and different MI calculation methods (I have tried 2 for tort and 5 for gaussian copula). Is this an issue with the pac_signals_tort function, or an issue within the permutation generation?

Hello again :-)

Just a quick question, in tensorpac.methods.ndpac, it seems that the PAC is squared before dividing by the number of points:

My understanding from reading the paper is that it should rather be something like:

pac = np.abs(np.einsum('i...j, k...j->ik...', amp, np.exp(1j * pha)))
s = pac**2
pac /= npts
# Set to zero non-significant values:
xlim = npts * erfinv(1 - p)**2
pac[s <= 2 * xlim] = 0.

Does that sound correct or did I miss something in the paper?

Thanks!