.. only:: html .. note:: :class: sphx-glr-download-link-note Click :ref:`here ` to download the full example code .. rst-class:: sphx-glr-example-title .. _sphx_glr_auto_examples_plot_spurious_pac.py: Spurious PAC with spikes ------------------------ This example demonstrates how a spike train can generate a significant amount of phase-amplitude coupling (PAC), even though there is no nested oscillations. This phenomenon is usually referred as spurious PAC. References ========== Dupre la Tour et al. (2017). Non-linear Auto-Regressive Models for Cross-Frequency Coupling in Neural Time Series. bioRxiv, 159731. .. code-block:: default import numpy as np import matplotlib.pyplot as plt import scipy.signal from pactools import Comodulogram, REFERENCES from pactools.utils.pink_noise import pink_noise from pactools.utils.validation import check_random_state def simulate_spurious_pac(n_points, fs, spike_amp=1.5, spike_fwhm=0.01, spike_fq=10., spike_interval_jitter=0.2, random_state=None): """Simulate some spurious phase-amplitude coupling (PAC) with spikes References ---------- Gerber, E. M., Sadeh, B., Ward, A., Knight, R. T., & Deouell, L. Y. (2016). Non-sinusoidal activity can produce cross-frequency coupling in cortical signals in the absence of functional interaction between neural sources. PloS one, 11(12), e0167351 """ n_points = int(n_points) fs = float(fs) rng = check_random_state(random_state) # draw the position of the spikes interval_min = 1. / float(spike_fq) * (1. - spike_interval_jitter) interval_max = 1. / float(spike_fq) * (1. + spike_interval_jitter) n_spikes_max = np.int(n_points / fs / interval_min) spike_intervals = rng.uniform(low=interval_min, high=interval_max, size=n_spikes_max) spike_positions = np.cumsum(np.int_(spike_intervals * fs)) spike_positions = spike_positions[spike_positions < n_points] # build the spike time series, using a convolution spikes = np.zeros(n_points) spikes[spike_positions] = spike_amp # full width at half maximum to standard deviation convertion spike_std = spike_fwhm / (2 * np.sqrt(2 * np.log(2))) spike_shape = scipy.signal.gaussian(M=np.int(spike_std * fs * 10), std=spike_std * fs) spikes = scipy.signal.fftconvolve(spikes, spike_shape, mode='same') noise = pink_noise(n_points, slope=1., random_state=random_state) return spikes + noise, spikes Let's generate the spurious PAC signal .. code-block:: default fs = 1000. # Hz n_points = 60000 spike_amp = 1.5 random_state = 0 # generate the signal signal, spikes = simulate_spurious_pac( n_points=n_points, fs=fs, random_state=random_state, spike_amp=spike_amp) # plot the signal and the spikes n_points_plot = np.int(1. * fs) time = np.arange(n_points_plot) / fs fig, axs = plt.subplots(nrows=2, figsize=(10, 6), sharex=True) axs[0].plot(time, signal[:n_points_plot], color='C0') axs[0].set(title='spikes + pink noise') axs[1].plot(time, spikes[:n_points_plot], color='C1') axs[1].set(xlabel='Time (sec)', title='spikes') plt.show() .. image:: /auto_examples/images/sphx_glr_plot_spurious_pac_001.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-script-out Out: .. code-block:: none /home/tom/work/github/pactools/examples/plot_spurious_pac.py:82: UserWarning: Matplotlib is currently using agg, which is a non-GUI backend, so cannot show the figure. plt.show() Compute the comodulograms .. code-block:: default methods = ['ozkurt', 'penny', 'tort', 'duprelatour'] low_fq_range = np.linspace(1., 20., 60) high_fq_range = np.linspace(low_fq_range[-1], 100., 60) low_fq_width = 4. # Hz # A good rule of thumb is n_surrogates = 10 / p_value. Example: 10 / 0.05 = 200 # Here we use 10 to be fast n_surrogates = 10 p_value = 0.05 n_jobs = 11 # prepare the plot axes n_lines = 2 n_columns = int(np.ceil(len(methods) / float(n_lines))) figsize = (4 * n_columns, 3 * n_lines) fig, axs = plt.subplots(nrows=n_lines, ncols=n_columns, figsize=figsize) axs = axs.ravel() for ax, method in zip(axs, methods): # compute the comodulograms estimator = Comodulogram(fs=fs, low_fq_range=low_fq_range, high_fq_range=high_fq_range, low_fq_width=low_fq_width, method=method, n_surrogates=n_surrogates, n_jobs=n_jobs) estimator.fit(signal) # plot the comodulogram with contour levels estimator.plot(contour_method='comod_max', contour_level=p_value, axs=[ax], titles=[REFERENCES[method]]) fig.tight_layout() plt.show() .. image:: /auto_examples/images/sphx_glr_plot_spurious_pac_002.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-script-out Out: .. code-block:: none [ ] 0% | 0.00 sec | comodulogram: ozkurt [ ] 2% | 7.11 sec | comodulogram: ozkurt [. ] 3% | 7.86 sec | comodulogram: ozkurt [.. ] 5% | 8.70 sec | comodulogram: ozkurt [.. ] 7% | 9.36 sec | comodulogram: ozkurt [... ] 8% | 10.14 sec | comodulogram: ozkurt [.... ] 10% | 10.90 sec | comodulogram: ozkurt [.... ] 12% | 11.58 sec | comodulogram: ozkurt [..... ] 13% | 12.24 sec | comodulogram: ozkurt [...... ] 15% | 12.78 sec | comodulogram: ozkurt [...... ] 17% | 13.39 sec | comodulogram: ozkurt [....... ] 18% | 13.89 sec | comodulogram: ozkurt [........ ] 20% | 14.51 sec | comodulogram: ozkurt [........ ] 22% | 15.04 sec | comodulogram: ozkurt [......... ] 23% | 15.76 sec | comodulogram: ozkurt [.......... ] 25% | 16.50 sec | comodulogram: ozkurt [.......... ] 27% | 17.15 sec | comodulogram: ozkurt [........... ] 28% | 17.78 sec | comodulogram: ozkurt [............ ] 30% | 18.27 sec | comodulogram: ozkurt [............ ] 32% | 18.89 sec | comodulogram: ozkurt [............. ] 33% | 19.82 sec | comodulogram: ozkurt [.............. ] 35% | 20.69 sec | comodulogram: ozkurt [.............. ] 37% | 21.52 sec | comodulogram: ozkurt [............... ] 38% | 22.49 sec | comodulogram: ozkurt [................ ] 40% | 23.41 sec | comodulogram: ozkurt [................ ] 42% | 23.98 sec | comodulogram: ozkurt [................. ] 43% | 24.44 sec | comodulogram: ozkurt [.................. ] 45% | 25.17 sec | comodulogram: ozkurt [.................. ] 47% | 26.03 sec | comodulogram: ozkurt [................... ] 48% | 26.59 sec | comodulogram: ozkurt [.................... ] 50% | 27.45 sec | comodulogram: ozkurt [.................... ] 52% | 28.73 sec | comodulogram: ozkurt [..................... ] 53% | 29.47 sec | comodulogram: ozkurt [...................... ] 55% | 30.71 sec | comodulogram: ozkurt [...................... ] 57% | 31.87 sec | comodulogram: ozkurt [....................... ] 58% | 32.54 sec | comodulogram: ozkurt [........................ ] 60% | 33.33 sec | comodulogram: ozkurt [........................ ] 62% | 33.69 sec | comodulogram: ozkurt [......................... ] 63% | 34.09 sec | comodulogram: ozkurt [.......................... ] 65% | 34.51 sec | comodulogram: ozkurt [.......................... ] 67% | 34.85 sec | comodulogram: ozkurt [........................... ] 68% | 35.20 sec | comodulogram: ozkurt [............................ ] 70% | 35.58 sec | comodulogram: ozkurt [............................ ] 72% | 35.90 sec | comodulogram: ozkurt [............................. ] 73% | 36.26 sec | comodulogram: ozkurt [.............................. ] 75% | 36.60 sec | comodulogram: ozkurt [.............................. ] 77% | 37.20 sec | comodulogram: ozkurt [............................... ] 78% | 37.73 sec | comodulogram: ozkurt [................................ ] 80% | 38.67 sec | comodulogram: ozkurt [................................ ] 82% | 39.54 sec | comodulogram: ozkurt [................................. ] 83% | 40.14 sec | comodulogram: ozkurt [.................................. ] 85% | 40.86 sec | comodulogram: ozkurt [.................................. ] 87% | 41.35 sec | comodulogram: ozkurt [................................... ] 88% | 41.92 sec | comodulogram: ozkurt [.................................... ] 90% | 42.73 sec | comodulogram: ozkurt [.................................... ] 92% | 43.46 sec | comodulogram: ozkurt [..................................... ] 93% | 44.04 sec | comodulogram: ozkurt [...................................... ] 95% | 44.66 sec | comodulogram: ozkurt [...................................... ] 97% | 45.22 sec | comodulogram: ozkurt [....................................... ] 98% | 46.41 sec | comodulogram: ozkurt [........................................] 100% | 47.01 sec | comodulogram: ozkurt [ ] 0% | 0.00 sec | comodulogram: penny [ ] 2% | 2.27 sec | comodulogram: penny [. ] 3% | 3.12 sec | comodulogram: penny [.. ] 5% | 3.77 sec | comodulogram: penny [.. ] 7% | 4.65 sec | comodulogram: penny [... ] 8% | 5.58 sec | comodulogram: penny [.... ] 10% | 6.36 sec | comodulogram: penny [.... ] 12% | 7.07 sec | comodulogram: penny [..... ] 13% | 8.12 sec | comodulogram: penny [...... ] 15% | 9.21 sec | comodulogram: penny [...... ] 17% | 10.47 sec | comodulogram: penny [....... ] 18% | 11.77 sec | comodulogram: penny [........ ] 20% | 12.68 sec | comodulogram: penny [........ ] 22% | 13.48 sec | comodulogram: penny [......... ] 23% | 14.27 sec | comodulogram: penny [.......... ] 25% | 15.02 sec | comodulogram: penny [.......... ] 27% | 15.74 sec | comodulogram: penny [........... ] 28% | 16.53 sec | comodulogram: penny [............ ] 30% | 17.21 sec | comodulogram: penny [............ ] 32% | 18.01 sec | comodulogram: penny [............. ] 33% | 18.81 sec | comodulogram: penny [.............. ] 35% | 19.61 sec | comodulogram: penny [.............. ] 37% | 20.45 sec | comodulogram: penny [............... ] 38% | 21.26 sec | comodulogram: penny [................ ] 40% | 22.14 sec | comodulogram: penny [................ ] 42% | 22.95 sec | comodulogram: penny [................. ] 43% | 24.00 sec | comodulogram: penny [.................. ] 45% | 24.77 sec | comodulogram: penny [.................. ] 47% | 25.49 sec | comodulogram: penny [................... ] 48% | 26.24 sec | comodulogram: penny [.................... ] 50% | 27.00 sec | comodulogram: penny [.................... ] 52% | 27.76 sec | comodulogram: penny [..................... ] 53% | 28.60 sec | comodulogram: penny [...................... ] 55% | 29.37 sec | comodulogram: penny [...................... ] 57% | 30.22 sec | comodulogram: penny [....................... ] 58% | 30.99 sec | comodulogram: penny [........................ ] 60% | 31.69 sec | comodulogram: penny [........................ ] 62% | 32.40 sec | comodulogram: penny [......................... ] 63% | 33.17 sec | comodulogram: penny [.......................... ] 65% | 33.87 sec | comodulogram: penny [.......................... ] 67% | 34.55 sec | comodulogram: penny [........................... ] 68% | 35.33 sec | comodulogram: penny [............................ ] 70% | 36.05 sec | comodulogram: penny [............................ ] 72% | 36.86 sec | comodulogram: penny [............................. ] 73% | 37.58 sec | comodulogram: penny [.............................. ] 75% | 38.35 sec | comodulogram: penny [.............................. ] 77% | 39.19 sec | comodulogram: penny [............................... ] 78% | 40.04 sec | comodulogram: penny [................................ ] 80% | 40.87 sec | comodulogram: penny [................................ ] 82% | 41.78 sec | comodulogram: penny [................................. ] 83% | 42.57 sec | comodulogram: penny [.................................. ] 85% | 43.97 sec | comodulogram: penny [.................................. ] 87% | 45.37 sec | comodulogram: penny [................................... ] 88% | 46.60 sec | comodulogram: penny [.................................... ] 90% | 47.84 sec | comodulogram: penny [.................................... ] 92% | 49.20 sec | comodulogram: penny [..................................... ] 93% | 50.54 sec | comodulogram: penny [...................................... ] 95% | 51.87 sec | comodulogram: penny [...................................... ] 97% | 53.21 sec | comodulogram: penny [....................................... ] 98% | 54.60 sec | comodulogram: penny [........................................] 100% | 55.68 sec | comodulogram: penny [ ] 0% | 0.00 sec | comodulogram: tort [ ] 2% | 2.58 sec | comodulogram: tort [. ] 3% | 3.58 sec | comodulogram: tort [.. ] 5% | 4.67 sec | comodulogram: tort [.. ] 7% | 6.27 sec | comodulogram: tort [... ] 8% | 8.03 sec | comodulogram: tort [.... ] 10% | 9.90 sec | comodulogram: tort [.... ] 12% | 11.70 sec | comodulogram: tort [..... ] 13% | 13.50 sec | comodulogram: tort [...... ] 15% | 15.24 sec | comodulogram: tort [...... ] 17% | 16.95 sec | comodulogram: tort [....... ] 18% | 18.92 sec | comodulogram: tort [........ ] 20% | 20.51 sec | comodulogram: tort [........ ] 22% | 22.48 sec | comodulogram: tort [......... ] 23% | 23.81 sec | comodulogram: tort [.......... ] 25% | 25.24 sec | comodulogram: tort [.......... ] 27% | 26.52 sec | comodulogram: tort [........... ] 28% | 28.08 sec | comodulogram: tort [............ ] 30% | 29.68 sec | comodulogram: tort [............ ] 32% | 31.39 sec | comodulogram: tort [............. ] 33% | 33.40 sec | comodulogram: tort [.............. ] 35% | 35.24 sec | comodulogram: tort [.............. ] 37% | 37.03 sec | comodulogram: tort [............... ] 38% | 38.90 sec | comodulogram: tort [................ ] 40% | 40.87 sec | comodulogram: tort [................ ] 42% | 42.95 sec | comodulogram: tort [................. ] 43% | 45.10 sec | comodulogram: tort [.................. ] 45% | 47.12 sec | comodulogram: tort [.................. ] 47% | 48.99 sec | comodulogram: tort [................... ] 48% | 51.03 sec | comodulogram: tort [.................... ] 50% | 52.87 sec | comodulogram: tort [.................... ] 52% | 54.69 sec | comodulogram: tort [..................... ] 53% | 56.80 sec | comodulogram: tort [...................... ] 55% | 58.72 sec | comodulogram: tort [...................... ] 57% | 60.65 sec | comodulogram: tort [....................... ] 58% | 62.69 sec | comodulogram: tort [........................ ] 60% | 64.75 sec | comodulogram: tort [........................ ] 62% | 66.72 sec | comodulogram: tort [......................... ] 63% | 68.68 sec | comodulogram: tort [.......................... ] 65% | 70.63 sec | comodulogram: tort [.......................... ] 67% | 72.58 sec | comodulogram: tort [........................... ] 68% | 74.84 sec | comodulogram: tort [............................ ] 70% | 76.63 sec | comodulogram: tort [............................ ] 72% | 78.80 sec | comodulogram: tort [............................. ] 73% | 81.34 sec | comodulogram: tort [.............................. ] 75% | 83.12 sec | comodulogram: tort [.............................. ] 77% | 84.51 sec | comodulogram: tort [............................... ] 78% | 86.28 sec | comodulogram: tort [................................ ] 80% | 87.57 sec | comodulogram: tort [................................ ] 82% | 88.95 sec | comodulogram: tort [................................. ] 83% | 90.20 sec | comodulogram: tort [.................................. ] 85% | 91.46 sec | comodulogram: tort [.................................. ] 87% | 93.09 sec | comodulogram: tort [................................... ] 88% | 94.34 sec | comodulogram: tort [.................................... ] 90% | 95.77 sec | comodulogram: tort [.................................... ] 92% | 97.54 sec | comodulogram: tort [..................................... ] 93% | 98.80 sec | comodulogram: tort [...................................... ] 95% | 100.35 sec | comodulogram: tort [...................................... ] 97% | 101.63 sec | comodulogram: tort [....................................... ] 98% | 103.23 sec | comodulogram: tort [........................................] 100% | 104.74 sec | comodulogram: tort [ ] 0% | 0.00 sec | comodulogram: DAR(10, 1) [ ] 2% | 0.07 sec | comodulogram: DAR(10, 1) [. ] 3% | 0.09 sec | comodulogram: DAR(10, 1) [.. ] 5% | 0.10 sec | comodulogram: DAR(10, 1) [.. ] 7% | 0.11 sec | comodulogram: DAR(10, 1) [... ] 8% | 0.13 sec | comodulogram: DAR(10, 1) [.... ] 10% | 0.14 sec | comodulogram: DAR(10, 1) [.... ] 12% | 0.16 sec | comodulogram: DAR(10, 1) [..... ] 13% | 0.18 sec | comodulogram: DAR(10, 1) [...... ] 15% | 0.19 sec | comodulogram: DAR(10, 1) [...... ] 17% | 0.21 sec | comodulogram: DAR(10, 1) [....... ] 18% | 0.23 sec | comodulogram: DAR(10, 1) [........ ] 20% | 0.28 sec | comodulogram: DAR(10, 1) [........ ] 22% | 0.30 sec | comodulogram: DAR(10, 1) [......... ] 23% | 0.47 sec | comodulogram: DAR(10, 1) [.......... ] 25% | 0.58 sec | comodulogram: DAR(10, 1) [.......... ] 27% | 0.66 sec | comodulogram: DAR(10, 1) [........... ] 28% | 0.76 sec | comodulogram: DAR(10, 1) [............ ] 30% | 0.84 sec | comodulogram: DAR(10, 1) [............ ] 32% | 0.92 sec | comodulogram: DAR(10, 1) [............. ] 33% | 1.00 sec | comodulogram: DAR(10, 1) [.............. ] 35% | 1.11 sec | comodulogram: DAR(10, 1) [.............. ] 37% | 5.54 sec | comodulogram: DAR(10, 1) [............... ] 38% | 5.61 sec | comodulogram: DAR(10, 1) [................ ] 40% | 5.75 sec | comodulogram: DAR(10, 1) [................ ] 42% | 5.88 sec | comodulogram: DAR(10, 1) [................. ] 43% | 6.03 sec | comodulogram: DAR(10, 1) [.................. ] 45% | 6.17 sec | comodulogram: DAR(10, 1) [.................. ] 47% | 6.34 sec | comodulogram: DAR(10, 1) [................... ] 48% | 6.50 sec | comodulogram: DAR(10, 1) [.................... ] 50% | 6.63 sec | comodulogram: DAR(10, 1) [.................... ] 52% | 6.77 sec | comodulogram: DAR(10, 1) [..................... ] 53% | 6.92 sec | comodulogram: DAR(10, 1) [...................... ] 55% | 11.76 sec | comodulogram: DAR(10, 1) [...................... ] 57% | 11.83 sec | comodulogram: DAR(10, 1) [....................... ] 58% | 11.97 sec | comodulogram: DAR(10, 1) [........................ ] 60% | 12.12 sec | comodulogram: DAR(10, 1) [........................ ] 62% | 12.25 sec | comodulogram: DAR(10, 1) [......................... ] 63% | 12.39 sec | comodulogram: DAR(10, 1) [.......................... ] 65% | 12.53 sec | comodulogram: DAR(10, 1) [.......................... ] 67% | 12.70 sec | comodulogram: DAR(10, 1) [........................... ] 68% | 12.83 sec | comodulogram: DAR(10, 1) [............................ ] 70% | 12.90 sec | comodulogram: DAR(10, 1) [............................ ] 72% | 12.94 sec | comodulogram: DAR(10, 1) [............................. ] 73% | 15.79 sec | comodulogram: DAR(10, 1) [.............................. ] 75% | 15.83 sec | comodulogram: DAR(10, 1) [.............................. ] 77% | 15.90 sec | comodulogram: DAR(10, 1) [............................... ] 78% | 15.96 sec | comodulogram: DAR(10, 1) [................................ ] 80% | 16.03 sec | comodulogram: DAR(10, 1) [................................ ] 82% | 16.13 sec | comodulogram: DAR(10, 1) [................................. ] 83% | 16.20 sec | comodulogram: DAR(10, 1) [.................................. ] 85% | 16.26 sec | comodulogram: DAR(10, 1) [.................................. ] 87% | 16.30 sec | comodulogram: DAR(10, 1) [................................... ] 88% | 16.38 sec | comodulogram: DAR(10, 1) [.................................... ] 90% | 16.46 sec | comodulogram: DAR(10, 1) [.................................... ] 92% | 19.70 sec | comodulogram: DAR(10, 1) [..................................... ] 93% | 19.72 sec | comodulogram: DAR(10, 1) [...................................... ] 95% | 19.78 sec | comodulogram: DAR(10, 1) [...................................... ] 97% | 19.81 sec | comodulogram: DAR(10, 1) [....................................... ] 98% | 19.88 sec | comodulogram: DAR(10, 1) [........................................] 100% | 19.95 sec | comodulogram: DAR(10, 1) [........................................] 100% | 25.40 sec | comodulogram: DAR(10, 1) /home/tom/work/github/pactools/examples/plot_spurious_pac.py:118: UserWarning: Matplotlib is currently using agg, which is a non-GUI backend, so cannot show the figure. plt.show() .. rst-class:: sphx-glr-timing **Total running time of the script:** ( 3 minutes 54.364 seconds) .. _sphx_glr_download_auto_examples_plot_spurious_pac.py: .. only :: html .. container:: sphx-glr-footer :class: sphx-glr-footer-example .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: plot_spurious_pac.py ` .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: plot_spurious_pac.ipynb ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_