3.1.1. Coupling Strength

In order to test the ability of our phase TE estimation method to detect phase-based directed interactions of varying intensity, we modify the coupling strength between the simulated signals, **x** and **y**, by varying the parameter *C*<sup>12</sup> in the range {0, 0.2, 0.5, 0.8}, with 0 indicating the absence of coupling and 0.8 a strong interaction between the two signals.

## 3.1.2. Noise and Signal Mixing

To asses the robustness of our proposal to realistic levels of noise and signal mixing, we do the following: we generate a noise time series *η*, with the same power spectrum of **x**, through the methodology proposed in [8]. Then, we add **x** and *η* to generate a noisy version of **<sup>x</sup>**, **<sup>x</sup>***<sup>η</sup>* = **<sup>x</sup>** + <sup>10</sup><sup>−</sup> SNR <sup>20</sup> *η*, where SNR is the signal to noise ratio. Likewise for **y**. Then, we mix **x***<sup>η</sup>* and **y***<sup>η</sup>* to simulate one of the effects of volume conduction, by doing **x***w <sup>η</sup>* = <sup>1</sup> <sup>−</sup> *<sup>w</sup>* 2 **x***<sup>η</sup>* + *<sup>w</sup>* 2 **y***η*, and **y***<sup>w</sup> <sup>η</sup>* <sup>=</sup> <sup>1</sup> <sup>−</sup> *<sup>w</sup>* 2 **<sup>y</sup>***<sup>η</sup>* <sup>+</sup> *<sup>w</sup>* 2 **x***η*, with *w* the mixing strength. We set the parameters SNR and *w* to 3 and 0.25 respectively, based on the results obtained in [8] for realistic values of noise and mixing for EEG signals. The coupling coefficient *C*<sup>12</sup> is held constant at a value of 0.5 to simulate couplings of medium strength.

#### 3.1.3. Narrowband Bidirectional Interactions

In this experiment, we aim to evaluate how our proposal deals with bidirectional interactions of localized frequency content. Particularly, we want to assess its performance for signals **x** and **y** containing a directed interaction from **x** to **y** at 10 Hz and an interaction in the opposite direction, from **y** to **x**, at 40 Hz. To generate such signals, first, we modify the model parameters of Area 2 so that it produces a signal **y** with a power spectrum peaking in the *γ* band [39]. The power spectrum of **x** remains as before. The coupling coefficient *C*<sup>12</sup> is again held constant at a value of 0.5. The change in the parameters of Area 2 leads to strong directed interactions from **x** to **y** around 10 Hz and 40 Hz. Then, we use a Morlet wavelet (Equation (10)) to filter both **x** and **y** at those frequencies (10 Hz and 40 Hz). The obtained real-valued narrowband time series are then combined as follows: **<sup>x</sup>**<sup>∗</sup> <sup>=</sup> **<sup>x</sup>**<sup>10</sup> *Hz* <sup>+</sup> **<sup>y</sup>**<sup>40</sup> *Hz* and **<sup>y</sup>**<sup>∗</sup> <sup>=</sup> **<sup>y</sup>**<sup>10</sup> *Hz* <sup>+</sup> **<sup>x</sup>**<sup>40</sup> *Hz*. Next, **<sup>x</sup>**<sup>∗</sup> and **<sup>y</sup>**<sup>∗</sup> are added to broadband noise generated following the same approach described in Section 3.1.2, with an SNR of 6.
