Phase-amplitude coupling

Authors: Soheila Samiee, Thomas Donoghue, Francois Tadel, Sylvain Baillet

This tutorial introduces the concept of phase-amplitude coupling (PAC) and the metrics used in Brainstorm to estimate it. These tools are illustrated first on simulated recordings. In separate tutorials, we illustrate how to use them on resting-state MEG recordings.

Contents

  1. Introduction
  2. Simulate signals
  3. PAC estimation with MVL
  4. Time-resolved PAC estimation with tPAC
  5. Further analysis of PAC files
  6. Practical suggestions for PAC analysis

Introduction

Phase-amplitude coupling (PAC)

The oscillatory activity in multiple frequency bands is observed in different levels of organization from micro-scale to meso-scale and macro-scale. Studies have been shown that brain functions are achieved with simultaneous oscillations in different frequency bands [Schutter and Knyazev, 2012]. Classical studies in this field were only focused on rhythms in each of these frequency bands, and it has been reported that these rhythms are linked to perception and cognition [Cohen, 2008]. However, it is revealed that not only examining brain activity in each single frequency band, but also the relation and interaction between oscillations in different bands, can be informative in understanding brain function. Thus, this concept increasingly received interest especially in the field of cognitive neuroscience. This interaction between several oscillations is also known as cross-frequency coupling (CFC).

Two forms of recognized CFC in brain rhythms are: phase amplitude coupling (PAC), and phase-phase coupling (PPC). In the first type, which is also called nested oscillations, the phase of the lower frequency oscillation (nesting) drives the power of the coupled higher frequency oscillation (nested), that results in synchronization of amplitude envelope of faster rhythms with the phase of slower rhythms. The second form is amplitude independent phase locking between n cycles of high frequency oscillation and m cycles of low frequency one. That's why it is also called n:m phase synchrony [Palva et al., 2005].

Among these two types, phase-amplitude coupling received more interests. It has been shown that behavioral tasks can modulate the phase amplitude coupling [Voytek et al., 2010], and also it is potentially involved in sensory integration, memory process, and attentional selection [Lisman and Idiart, 1995, Lisman, 2005, Schroeder and Lakatos, 2009]. This coupling is observed in several brain regions including hippocampus, basal ganglia, and neocortex; and these observations are reported in rats, mice, sheep, and monkeys, as well as humans [Tort et al., 2010].

Following figure shows a schematic of phase amplitude coupling. In the top signal, we have the sum of a fast and slow oscillations, where the power of fast oscillation's envelope changes with the phase of the slower oscillation, which is a sample of PAC. The bottom signal shows only the fast oscillation and the variation in its power. As it is obvious from comparison of two signals, the fast rhythm's power is always maximum, at a certain phase of slower oscillation, this phase is called coupling phase.

pac_schematic2.jpg

Measures of cross frequency phase amplitude coupling can monitor the relationship between the activities that modulate low frequency oscillations like sensory or motor inputs, and the local cortical activities such as local computations that are correlated to amplitude of higher frequency oscillation [Canolty and Knight, 2010].

All these interesting features of this coupling resulted in proposing several methods for its measuring. Each of these methods has certain limitations and also advantages over the others, and can be used for a particular purpose; that's why no preferred standard method has been chosen for this estimation yet [Tort et al., 2010]. One of these methods, which is implemented in Brainstorm, is called Mean Vector Length (MVL), proposed by Canolty et al. (2006).

On the other hand, since phase-amplitude coupling changes over time with interinsic events and external stimuli, time-resolved estimation of this phenomena received increasing interest. Among suggested methods for estimating the dynamics of phase-amplitude coupling, tPAC [Samiee & Baillet, 2017], is implemented in Brainstorm.

These methods (MVL and tPAC) and the step by step instruction of using them are described in this section of the tutorial.

Two measures of PAC:

1. Time-resolved PAC estimation: tPAC

This approach in principle searches for the $$f_P$$ oscillation with strongest PAC to $$f_A$$ bursts, over a time-window, which slides on the input electrophysiological data. Following figure, summarizes the steps of the algorithm (see [Samiee & Baillet, 2017] for more details).

Fig1_method_Final.jpg

2. Mean Vector Length (Modulation Index)

Canolty et al. (2006) pointed out that a time series defined in the complex plane by $$A_{f_A}. e^{i \phi_{f_p}}$$ could be used to extract a phase-amplitude coupling measure. In this formula $$A_{f_A}$$ is the envelope of fast oscillation, and $$\phi_{f_p}$$ is the phase of slow oscillation.

Therefore, after filtering in fast and slow oscillation, and extracting the phase of slow, and the amplitude of fast rhythm; each instantaneous fast oscillation amplitude component in time is represented by the length of the complex vector, whereas the slow oscillation phase of the same time point is represented by the vector angle (see following figure).

mvl_step3.jpg mvl_steps.jpg

At the absence of phase-amplitude coupling, the plot of the $$A_{f_A} .e^{i \phi_{f_p}}$$ time series in the complex plane is characterized by a roughly uniform circular density of vector points, symmetric around zero, because the $$A_{f_A}$$ values (averaged over cycles of slow oscillation) are approximately the same for all phases. If there is modulation of the $$f_A$$ amplitude by the $$f_P$$ phase, the $$A_{f_A}$$ would be higher at certain phases than others. This higher amplitude for certain angles will lead to a “bump” in the polar plot of the $$A_{f_A}.e^{i \phi_{f_p}}$$, leading to loss of symmetry around zero. This loss of symmetry can be inferred by measuring the length of the vector obtained from the mean over all points in the complex plane. It is thus assumed that a symmetric distribution as it occurs during lack of coupling leads to a small mean vector length (because the points in the different phases would cancel each other), whereas the existence of coupling leads to a larger mean vector length [Tort et al. 2010]. For more detail read [Canolty et al, 2006].

References

Canolty RT, Edwards E, Dalal SS, Soltani M, Nagarajan SS, Kirsch HE, Berger MS, Barbaro NM, Knight RT (2006). High gamma power is phase-locked to theta oscillations in human neocortex, Science, 313(5793), 1626-1628.

Canolty RT, Knight RT (2010). The functional role of cross-frequency coupling. Trends in cognitive sciences, 14(11), 506-515

Cohen MX (2008). Assessing transient cross-frequency coupling in EEG data. Journal of neuroscience methods, 168(2), 494-499.

Palva JM, Palva S, Kaila K (2005). Phase synchrony among neuronal oscillations in the human cortex. The Journal of Neuroscience, 25(15), 3962-3972.

Samiee, Soheila, and Sylvain Baillet (2017). Time-resolved phase-amplitude coupling in neural oscillations. NeuroImage, 59, 270-279.Schutter DJ, Knyazev GG(2012). Cross-frequency coupling of brain oscillations in studying motivation and emotion. Motivation and emotion, 36(1), 46-54.

Tort AB, Komorowski R, Eichenbaum H, Kopell N (2010). Measuring phase-amplitude coupling between neuronal oscillations of different frequencies. Journal of neurophysiology, 104(2), 1195-1210.

Voytek B, Canolty RT, Shestyuk A, Crone NE, Parvizi J, Knight RT (2010). Shifts in gamma phase–amplitude coupling frequency from theta to alpha over posterior cortex during visual tasks.Frontiers in human neuroscience, 4.

Simulate signals

For PAC analysis you can use your own data, or generate a synthesized data containing cross-frequency phase-amplitude coupling, with your preferred parameters. Here we explain how to produce simulated data. The model used for data generation is a simple method introduced in [Tort et al. 2010]. In this section, we generate a dataset of synthesized signal and analyze it with available phase-amplitude coupling estimation tool in brainstorm.

Process options

PAC estimation with MVL

Input data

You can use your own data set, or the signal that is generated with brainstorm simulator. Here we used the simulated signal that is generated with the steps explained in the previous section.

PAC estimation

To extract the phase-amplitude coupling from this signal:

Process options

Input options:

Estimator options:

Loop options: Should be left at default options unless you know how to use them.

Output options:

File contents

The files saved by this process have the same structure as the time-frequency files, with an additional "sPAC" field. To review its contents, right-click on the PAC file > File > View file contents.

Some of the relevent fields in the maxPAC files:

Recommendations

PAC estimation using the MVL algorithm

In order to have a reliable result from this method it is required to use a signal which length is at least ten cycles of the slowest oscillation in your low oscillation band.

Considering this point is more important in analysis of real databases, where the noise level (and/or background brain activity) can be higher than synthesized data, and the coupling intensity can be low. Thus, if you want to examine the coupling for slow oscillations in [2, 14] Hz, it would be better to use a signal with minimum length of 10 cycles of the slowest oscillation, which would be 10 x 0.5 = 5 S.

Time-resolved PAC estimation with tPAC

Input data

For tPAC estimation you can use your own data set or simulated data. For illustration here, we simulated three time series with phase-amplitude coupling using brainstorm tool and combine them. The script to reproduce such signal can be found here. The resulting time series looks like:

simulated_data.png

During the first half of the signal (10 s duration), the phase of slow oscillation at 9 Hz is coupled to the amplitude of a faster rhythm at 115 Hz. In the second half, the first coupling mode is terminated and the average of two other simultaneous modes appears. These modes are $$f_{P_2}$$ = 13 Hz, $$f_{A_2}$$ = 145 Hz, and $$f_{P_3}$$=5 Hz, $$f_{A_3}$$ =87 Hz, respectively. The signal-to-noise ratio is set to 6 dB, and the preferred coupling phase in the three modes are 270, 0, and 180°, respectively. The coupling strength is kept constance and identical in all modes.

Time-resolved PAC

To extract the phase-amplitude coupling from this signal:

Process options

Input options:

Estimator options:

Output options:

File contents

The files saved by this process have the same structure as the time-frequency files, with an additional "sPAC" field. To review its contents, right-click on the PAC file > File > View file contents.

tPAC_output_file.png

Note: Preferred phase of coupling extracted for each pair of $$f_P$$-$$ f_A$$ is saved under sPAC.DynamicPhase. It can be exported to matlab and visualized there.

Comodulogram of time-resolved PAC:

To extract the comodulogram from the time-resolved PAC:

Process Options

Input options:

Method:

Time-resolved PAC in time - $$f_P$$ domain

...

Further analysis of PAC files

All functions for further analysis of PAC files (e.g., deriving the average of PAC maps, producing comodulograms, running a z-score transform with respect to baseline levels or surrogate data, etc. ) are available as processed located in Frequency > Time-resolved Phase-amplitude coupling.

PAC_further_analysis.png

Note: Basic functions for averaging or running statistics are not yet avalaible with Brainstorm [working on it...]. Yet, PAC output files are stored as Matlab files and contents can be readily processed with custom scripts of your own.

Yet, surrogate analysis is only available for single band tPAC maps, and can be accessed here: PAC_surrogate.png

Further, tPAC traces can be z-scored with respect to surrogate data using the Process2 tab.

Practical suggestions for PAC analysis

How to best get started with PAC analysis? How do you decide about setting parameters? Here, we provide basic recommendations that we hope can be helpful:

General suggestions

Time-resolved analysis

For time-resolved analysis, if you are not sure about the exact frequencies that you want to focus on, one practical pipeline is to extract comodulogram with a relatively wide band (based on your initial hypothesis), and then find the main mode of coupling from the comodulogram. In the next step, focus on that particular mode of coupling to evaluate dynamical changes in coupling intensity (or preferred phase of coupling) with a single band analysis (single band option does not provide comodulogram but reflect a trace of dynamical changes in coupling parameters). With this strategy the temporal resolution will be optimized for the main mode of coupling, and also the computational load for intensive analysis (and checking significance of outcomes) will be much lower than looking at the whole interval.

Following steps provide more details:

  1. Start with comodulogram and find the main mode (or modes) of coupling based on comodulogram.

  2. To decide about length of sliding time window in this step, you can start with two cycles of the slowest oscillation in the band of interest for $$f_P$$. For example, if you are interested in investigating potential PAC between theta to alpha oscillations (4-12 Hz) for $$f_P$$ and gamma rhythms for $$f_A$$, a conservative short sliding time window can be set to two cycles of 4 Hz (i.e. $$ 2 \times \frac{1}{4 Hz} = 2 \times 0.25 s = 0.5 s = 500 ms$$. Certainly, if you have long inputs increasing the length of the window would improve the SNR.

The only point that you should keep in mind is that, tPAC do sparse estimation of PAC in a 3D space of fP, fA and time; therefore, you need enough number of points in time dimension to be able to retrieve a reasonable projection of it on fp - fA plain (called comodulogram). Hence, the maximum length of the window can be selected based on the input length. For example, in the case where your data is around 5 minutes long, if you set the window length to 20 s, considering the 50% overlap of windows you will end up with 5 * 60s / 10s = 30 points along time dimension (which is used in averaging for comodulogram extraction). It is not recommended to go below a certain number (e.g 30).

  1. For extracting comodulogram you can also use none-time-resolved option (e.g. MVL) which use the whole time window and does not divide it into smaller chunks. It should be noted that in one experiment we observed higher resolution from comodulogram extracted using time-resolved method compared to a traditional option (see Fig. 5 in [Samiee & Baillet 2017]: https://www.sciencedirect.com/science/article/pii/S1053811917306195#fig5). The reason behind this improvement could be similar to the improvement in power spectral density estimation using Welch approach compared to original Bartlett's method.

  2. If you can find an obvious mode of coupling in the comodulogram, it would be better to investigate that mode of coupling using single band analysis for all regions of interest and time intevals (rather than maxPAC value extracted from comodulogram).
  3. If there is no obvious mode of coupling in averaged comodulogram, it is not recommended to go further in analysis. It would be better to take a more close look at modes availble in each single recording, and consistancy of preferred phase of coupling in time and in different trials/subjects.

  4. If you observe a blub in comodulogram which looks like a vertical line, you need to further investigate your data to make sure this observed coupling is not a pseudo coupling coming from sharp event artifacts (e.g. spikes in LFP) or asymmetric oscillations (e.g. oscillations with saw tooth waveform shape). See [Cole & Voytek 2017; Aru et al. 2015] for more details.

A quick summary of the most important practical points in PAC analysis

1. Temporal resolution and minimum $$f_P$$

2. Suggestion for deciding about initial $$f_P$$ band

3. Coupling strength for the main coupled pair




Feedback: Comments, bug reports, suggestions, questions
Email address (if you expect an answer):


Tutorials/TutPac (last edited 2023-04-19 19:00:56 by RaymundoCassani)