Global synchronization measure applied to brain signals
data

Xhilda Dhamo¹[0000-0002-5157-7075], Eglantina Kalluçi¹[0009-0001-9039-1310], Gérard Dray²[0000-0003-1525-5682], Coralie Reveille ³[0000-0002-2511-6611], Arnisa Sokoli¹[0000-0002-4039-2024], Stefan Perrey³[0000-0002-8741-629X], Gregoire Bosselut³[0000-0002-0395-0998], Stefan Janaqi²[0000-0002-6536-3799]

¹ Department of Applied Mathematics, Faculty of Natural Sciences, University of Tirana,
Tirana, Albania

² EuroMov Digital Health in Motion, University of Montpellier, IMT Mines Ales, Ales, France

³ EuroMov Digital Health in Motion, University of Montpellier, IMT Mines Ales, Montpellier,
France

xhilda.merkaj@fshn.edu.al

**Abstract.** We investigate a method to assess brain synchronization in individuals who fulfill a cooperation task. Our input is a couple of signals from functional Near-Infrared Spectroscopy Data Acquisition and Pre-processing technology that is used to capture the brain activity of an individual by measuring the oxyhemoglobin (HbO) level. Then, we use the visibility graph approach to map each HbO signal into a network. We estimate the signal synchronization by studying a global measure, related to eigenvalues of Laplacian matrix, in each constructed visibility graph. We consider the autonomous evolution of one isolated node to be a Rössler function. Then, the synchronization of signals can be characterized by a little number of parameters that could be employed to classify the sources of signal. Unlike prior research in this area, our aim is to examine the circumstances in which synchronization occurs in various individuals and within different hemispheres of the prefrontal cortexes of the same individual. Experimental results show that the conditions for synchronization vary in different individuals, and they are different even for the distinct prefrontal cortical hemispheres of the same individual.

**Keywords:** visibility graph, Laplacian, synchronization measure, brain synchronization

# 1 Introduction

Network Science comes out as an efficient tool for identifying, representing and predicting distinct collective phenomena in numerous complex systems (Newman, 2018, Halvin et al. 2012, Chen et al. 2014, Arenas et al. 2008, Boccaletti et al. 2006, Barrat et al. 2008). Synchronization problem appears to be a hottest topic in dynamic processes in networks (Ding et al. 2019). It is defined as a process wherein many systems adjust
a given property of their motion due to a suitable coupling configuration, or to an ex-
ternal forcing (Boccaletti et al. 2006). Synchronization processes are everywhere in
nature and they have been applied to a wide variety of areas such as biology, ecology,
climatology, sociology, technology, smart grid, secure communication, neurology and
many more (Pikovsky et al. 2001, Osipov et al. 2007, Arenas et al. 2008, Dörfler &
Francesco 2014, Jiruska et al. 2013).

Visibility graphs present a powerful technique to study time series in the context of
complex networks (Lacasa et al. 2008). These kinds of graphs were established as a
method to map time series into networks (Lacasa et al. 2008) with the purpose of mak-
ing use of the tools of network science to describe the structure and dynamics of time
series. Explorations on visibility graphs are mainly focused on two different directions:
(i) canonical dynamics such as stochastic and chaotic processes (Brú et al. 2017, Gon-
çalves et al. 2016, Lacasa et al. 2009, Luque et al 2013, Luque et al. 2011); (ii) a feature
extraction procedure to make statistical learning (Bhaduri & Ghosh 2016, Hou et al.
2016, Long et al. 2014). The implementation of visibility graph to neuroscience is in
its beginning and has been limited to the analysis of the electroencephalogram (EEG)
(Mira- Iglesias et al. 2016, Bhaduri & Ghosh 2016) and functional magnetic resonance
imaging (fMRI) data (Sannino et al. 2017).

In this paper we use visibility graphs to map functional Near-Infrared Spectroscopy Data Acquisition and Pre-processing (fNIRS) time series into networks and then study the synchronization dynamics in the constructed networks. fNIRS is a technology used to measure brain activity (Li et al. 2020). The HbO signals are captured by the optods positioned in left and right hemispheres of the prefrontal cortex (PFC) of the participants in the experiment. Then, we identify a few parameters that permit to characterize the synchronization process. These parameters could be used further to classify the sources of signal. The analysis of the synchronization phenomena in these signals show that the synchronization does not happen under the same conditions for different participants in the experiment, but even for different PFC hemispheres of one participant. However, these are preliminary results and need to be compared to alternative methods.

Recent research articles (Li et al. 2021, Wang et al. 2022) have outlined an approach to studying the synchronization of individuals’ brain when collaborate to complete a particular task. They propose using the Wavelet Transform Coherence (WTC) to identify locally phase locked behavior between two time-series by measuring cross-correlation between the time series as a function of frequency and time. Furthermore, they used sliding windows approach and k-means clustering to obtain a set of representative inter-brain network states during different communication tasks. Our experiment is performed only in one communication task and our goal in this paper is to compare the conditions under which synchronization is achieved in different individuals and also in different PFCs of the same participant.

The rest of the paper is organized as follows. In the second section we define pre-
liminaries on network theory, introduce the visibility graph and its properties and give
the mathematics behind the synchronization measure through undirected, unweighted
networks. The third section describes the generation of the data used in this study. In
addition, we provide the reader with results obtained when studying and analyzing synchronization dynamics in brain activity data. Conclusion summarizes once more all the work conducted and results obtained from our analysis.

# 2 Background and methods

## 2.1 Preliminaries in networks

Throughout this paper we will refer to a graph (network) as a pair $G = (V, E)$ where $V$ is called the vertex set and $E$ is the edge set. This study is focused only on undirected and unweighted networks with symmetric adjacency matrix $A$ whose entries are $a_{ij} = 1$ if there exist a link between nodes $i$ and $j$ and 0 otherwise. The degree of one node $i$ is defined as $k_i = \sum_{j=1}^{|V|} a_{ij}$, $i = 1, 2, ..., |V|$ and the Laplacian matrix associated to $G$ is $L = D - A$, where $D = \text{diag}\{k_1, k_2, ..., k_{|V|}\}$. It is known that $L$ is positive semi-definite matrix and all its eigenvalues are real and non-negative. In addition, they can be ordered as $0 = \lambda_1 < \lambda_2 < \lambda_3 < \cdots < \lambda_N$, where $\lambda_2$ is known as the algebraic connectivity of the network (Newman 2018, Estrada and Knight 2015, Barabási & Pósfai 2016, Boccaletti et al. 2006).

## 2.2 The visibility graph

The procedure to construct a visibility graph is described in detail in (Lacasa et al. 2008, 2009, 2012). Let's consider a time series with $N$ data measured at times $t_i, i = 1, 2, ..., N$ with values $x_i, i = 1, 2, ..., N$ and consecutive time points $(t_i, x_i)$, $(t_k, x_k)$ and $(t_j, x_j)$. Time points $(t_i, x_i)$ and $(t_j, x_j)$ are visible and consequently will become two connected nodes in the visibility graph if for any point $(t_k, x_k)$ between them, they fulfill the following inequality:

$$x_k < x_j + (x_i - x_j) \frac{t_j - t_k}{t_j - t_i} \quad (1)$$

The network constructed by the above condition has four main properties: it is connected, undirected, invariant under affine transformations of the series data and it can be applied to every kind of time series (Lacasa et al. 2008). The construction of visibility graph is illustrated schematically in figure 1 given a time series with $N=20$.
<figure><img src="image_1.png" /><figcaption>Fig. 1. Construction of visibility graph corresponding to a univariate time series. Adapted from Lacasa et al. 2008.</figcaption></figure>

## 2.3 Synchronization dynamic

Synchronization is the process where many systems adjust a given property of their motion due to a suitable coupling configuration, or to an external forcing (Arenas et al. 2008, Boccaletti et al. 2014, Barrat et al. 2008).

Let's consider a system (connected network) formed by N identical m dimensional dynamical system (oscillators) whose states are represented by the vector $X = \{x_1, x_2, ..., x_N\}$ where each of $x_i$ is a m-dimensional vector ($x_i$ stands for the nodes in the network). The evolution of the vector field $x_i$ could be determined by the ordinary differential equation $\dot{x}_i = f(x_i) (f: R^m \rightarrow R^m)$. The equation of motion is defined as:

$$ \dot{x}_i = f(x_i) - \sigma \sum_{j=1}^{N} L_{ij} h(x_j) \quad (2) $$

where $L$ is the Laplacian matrix, $h: R^m \rightarrow R^m$ is a vectorial output function and $\sigma$ is the coupling strength. The completely synchronized state of this network with identical dynamics is computed as the solution of equation (2) where $x_1(t) = x_2(t) = \dots = x_N(t) = s(t)$. In this synchronized state, $s(t)$ is also solution of the equation $\dot{s} = f(s)$. This subspace in the state space of equation (2) where all the oscillators evolve synchronously on the same solution of the isolated oscillator $f$ is called the synchronization manifold. Once, all the oscillators are set at the synchronization manifold, they will remain synchronized and the most important topic now is to evaluate the stability of the synchronized manifold in presence of small perturbations $\delta x_i$. The synchronized manifold is stable if the perturbations decay in time, otherwise it is not stable.

The stability of the synchronized manifold is measured making use of the master stability function. We consider the evolution of small $\delta x_i$ as linear and write $x_i = s(t) + \delta x_i$. Furthermore, we write the Taylor series for the functions $f$ and $h$ as $f(x_i) = f(s) + Jf(s)\delta x_i$ and $h(x_i) = h(s) + Jh(s)\delta x_i$, where $Jf(s)$ and $Jh(s)$ are the Jacobian matrices of the functions $f$ and $h$ on $s$ and obtain the variational equations for $\delta x_i$:
$$ \delta \dot{x}_i = Jf(s)\delta x_i - \sigma Jh(s) \sum_{j=1}^{N} L_{ij}\delta x_j \quad (3) $$

Then, $\delta x$ is projected into the eigenspace spanned by the eigenvectors of the Laplacian matrix and the linear variational equations are written below:

$$ \dot{\xi}_l = [Jf(s) - \sigma\lambda_l Jh(s)]\xi_l \quad (4) $$

Where $\xi_l$ is the eigenmode corresponding to the eigenvalue $\lambda_l$ of L. Since, all the variational equations have the same form, but only differ from the term $\alpha = \sigma\lambda_l$, we can write the variational equations in vectorial form:

$$ \dot{\xi} = [Jf(s) - \alpha Jh(s)]\xi \quad (5) $$

To evaluate the stability of the synchronization manifold, one computes the largest Lyapunov exponent $\lambda_{max}$ as a function of $\alpha$ based on the master variational equation. The synchronization manifold will be stable for all $\alpha$, which give a negative value for $\lambda_{max}$. This approach is given in detail by (Pecora and Carroll, 1998). Therein one finds a better understanding of the Master Stability Function (MSF) and its three possible classes. Here, our intention is to determine $\alpha_1$ and $\alpha_2$, so that $\sigma\lambda_l > \alpha$ ($\alpha_1 < \sigma\lambda_l < \alpha_2, l = 2, ..., N$) in case of class II (III) MSF. Since, the eigenvalues of the Laplacian matrix are non-negative we can obtain the following inequalities:

$$ \sigma\lambda_2 \le \sigma\lambda_3 \le \dots \le \sigma\lambda_N > \alpha \quad (\alpha_1 < \sigma\lambda_2 \le \sigma\lambda_3 \le \dots \le \sigma\lambda_N < \alpha_2) \quad (6) $$

The synchronization manifold s will be stable only for:

$$ \sigma > \frac{\alpha}{\lambda_2} \quad \left( \frac{\lambda_2}{\lambda_N} < \sigma < \frac{\alpha_2}{\alpha_1} \right) \quad (7) $$

The synchronization error is computed as:

$$ E = \lim_{T \to \infty} \int_{0}^{T} \sum_{j \neq 1} \|x_j(t) - x_1(t)\| \quad (8) $$

# 3 Experimental Results

## 3.1 Experiment setup

One cooperative task called “MapTask” is given to two participants, which together create a dyad. The guide participant (pA) had 35 icons and a path drawn on its screen. The follower participant (pB) had exactly the same icons but didn’t have the path drawn on its screen. The target of the experiment is that at the end of the task, pB should have drawn the same path as shown on the pA’s screen using the arrows in the keyboard and the instructions given by pA.
<figure><img src="image_2.png" /><figcaption>Fig. 2. Cooperative task “MapTask” and HbO signals of one dyad</figcaption></figure>

The technology used to measure brain activity is Functional Near-Infrared Spectroscopy Data Acquisition and Pre-processing (fNIRS). Each of the participants had two optods in their PFC: one positioned in the left hemisphere (hL) and the other in the right hemisphere (hR). These optods captured the oxyhemoglobin (HbO) and deoxyhemoglobin (HbR) signals. Considering that the HbO signal is more sensitive to changes in cerebral blood flow than the HbR signal, we focused on the HbO signal (Wang et al. 2022). 18 participants took part in the experiment and they were divided in 9 dyads. We consider the signal in the first five minutes of the beginning of the experiment. The cooperation task given to the participants in the experiment is illustrated in figure 2. A HbO signal is given too.

## 3.2 Brain synchronization

The HbO time series are mapped into visibility networks as described in section 2. All the network constructed have the same number of nodes (3148 nodes) which corresponds to the time points in the first five minutes of the duration of the experiment. There are two networks for each participant, corresponding to the signals obtained by measuring the HbO in left (IPFC) and right prefrontal cortex hemispheres (rPFC). To study the synchronization in each network we have considered as the autonomous evolution function of one isolated node the Rössler function $f=[-y-z;x+0.2y;0.2+z(x-9.0)]$ and as the output function $h = [0, y, 0]$, which determines a class II MSF. For simplicity, from now on we will refer as pAhL (pAhR) the IPFC (rPFC) HbO signals of pA and as pBhL (pBhR) the IPFC (rPFC) HbO signals of pB.

The results showed that different participants reach the synchronization state for different conditions. Furthermore, even the two PFC of the same participant do not reach synchronization under the same conditions. Here we represent each participant by the algebraic connectivity of its corresponding visibility graph. The values of the parameter decrease with the increasing values of the algebraic connectivity. The values of the parameter $\sigma$ for which the synchronization state is stable and the synchronization errors are reported in figure 3 for pAhL and pAhR signals and in figure 4 for pBhL and pBhR signals.
<figure><img src="image_3.png" /><figcaption>**Fig. 3.** A: The values of parameter $\sigma$ when MLE becomes negative for pAhL. Inset illustrates the same for pAhR. B: The synchronization error for all pAhL participants in the main figure and for pAhR in the inset. Main plot' ('Inset') legend corresponds to the main (insets) plots in A and B. Red dashed line is used here to identify when MLE becomes negative.</figcaption></figure>

<figure><img src="image_4.png" /><figcaption>**Fig. 4.** A: The values of parameter $\sigma$ when MLE becomes negative for pBhL. Inset illustrates the same for pBhR. B: The synchronization error for all pBhL participants in the main figure and for pBhR in the inset. 'Main plot' ('Inset') legend corresponds to the main (insets) plots in A and B. Red dashed line is used here to identify when MLE becomes negative.</figcaption></figure>

The error of synchronization (Ding et al. 2020) is computed for each oscillator as
$$E_i(t) = [x_i(t) - s_i^x(t), y_i(t) - s_i^y(t), z_i(t) - s_i^z(t)], \quad i = 1, \dots, N \text{ and the total}$$
synchronization error is computed as $E(t) = \sqrt{\frac{1}{N} \sum_{i=1}^{N} E_i^T E_i}$.

In figure 5, there are two participants of the first dyad and their IPFC and rPFC HbO signals.

<figure><img src="image_5.png" /><figcaption>**Fig. 5.** Time needed for one dyad to be synchronized considering both PFC hemispheres of each participant for all oscillators. **pAhL** refers to the IPFC of participant pA. The values of parameter $\sigma = 13.48$ (from now when MLE becomes negative) and $\lambda_2 = 0.0186$. **pAhR** refers to the rPFC of participant pA. The values of parameter $\sigma = 14.41$ and $\lambda_2 = 0.0173$. **pBhL** refers to the IPFC of participant pB with $\sigma = 127.03$ and $\lambda_2 = 0.002$. **pBhR** refers to the rPFC of participant pB with $\sigma = 111.68$ and $\lambda_2 = 0.0022$.</figcaption></figure>

Figure 5 indicates that the higher the value of algebraic connectivity, more strongly connected is the network and synchronization state is reached faster. In this dyad, the pA participant reaches synchronization faster than the pB participant. If we focus on one participant both PFC hemispheres need approximately the same time to reach synchronization. One different case is reported in figure 6, taking into consideration another dyad, when the two PFC hemispheres of one participant need different amounts of time to be synchronized.
<figure><img src="image_6.png" /><figcaption>**Fig. 6.** Time needed for both PFC hemispheres of one guider participant to be synchronized. **pAhL** refers to the lPFC of participant pA. The values of parameter σ = 172.4 and λ₂ = 0.0015. **pAhR** refers to the rPFC of participant pA. The values of parameter σ = 79.55 and λ₂ = 0.0031.</figcaption></figure>

The total synchronization error for all dyads is illustrated in figure 7.

<figure><img src="image_7.png" /><figcaption>**Fig. 7.** Total synchronization error. Dyad 1 is the same dyad presented in figure 5 and Dyad 2 refers to the same participants presented in figure 6. The legend is the same for all plots.</figcaption></figure>

We can identify three classes of HbO synchronization in the experiment conducted: (i) the HbO signals are approximately synchronized at the same time (dyad 6 and 9 despite the pBhR signal and dyad 3 despite the pAhL); (ii) the HbO signals of the pA are synchronized faster than the ones of pB (dyad 1); (iii) the HbO signals of pB are synchronized faster than the ones of pA (dyad 5, 7 and 8). Comparing the HbO synchronization of the hemispheres of the PFC it is noticed that in half of the participants the synchronization is reached at the same time for both hemispheres (dyad 1, 3, 5, 6, 8).
Furthermore, results indicate that if two individuals interacting as one dyad do not synchronize at approximately the same time, regardless of whether their two PFCs synchronize simultaneously, they require more time to complete the experiment (dyad 1 and dyad 5 completed the experiment in 10.48 and 22.1 minutes). In cases when the PFCs of different participants partnering within one dyad achieve synchronization at approximately the same time, the experiment is finished in a shorter time (dyad 2 and dyad 9 completed the experiment in 8.1 and 8.5 minutes). In other cases, when synchronization is achieved at approximately the same time only for two PFCs of different participants collaborating within a dyad, the duration of the experiment is larger than in the previous scenario (dyad 3, 4, 6, 7, 8 accomplished the experiment in 11.9, 14.8, 14.9, 17.9 and 12 minutes respectively).

## 4 Conclusions

In this paper we have investigated the brain synchronization problem from the perspective of complex networks. The human brain activity was measured using the fNIRS technology during the ‘MapTask’ experiment. There were 18 participants who took part in the experiments, grouped in 9 dyads. Two optods positioned in the left and right hemispheres of the prefrontal cortex were used to capture the HbO signals. Furthermore, we have used the visibility graph approach to convert the HbO time series into networks, where each node represents one time point and two nodes are linked if their corresponding time point has visibility with each other. In addition, we have used the constructed networks to study the global synchronization phenomena.

Experimental results indicated that different participants reached the synchronization state under different conditions. Furthermore, we have observed that even the two prefrontal cortexes of the same participant do not reach synchronization for the same conditions. The total synchronization error displayed in figure 7 suggested that in 8 participants (dyad 1, 3, 5, 6, 8) both prefrontal cortex hemispheres reach synchronization approximately at the same time. This work pointed out that when PFCs of different participants collaborating within same dyad synchronize roughly at the same time, the experiment is finished faster than in cases when they do not fulfil this property, even though both PFCs of same individual synchronize at the same time.

## Acknowledgement

We would like to thank Guilhem Belda, President of the Semaxone company (Rochefort-du-Gard, France), for his help in developing the application and providing the audio equipment.
References

1. Arenas, A., Díaz-Guilera, A., Kurths, J., Moreno, Y., Zhou, C.: Synchronization in complex networks. Physics Reports. 469, 93-153 (2008). https://doi.org/https://doi.org/10.1016/j.physrep.2008.09.002

2. Barabási, A.-L., Pósfai, M.: Network Science. Cambridge University Press (2016)

3. Barrat, A., Barthélemy, M., Vespignani, A.: Dynamical Processes on Complex Networks. Cambridge University Press (2008)

4. Bhaduri, A., Ghosh, D.: Quantitative Assessment of Heart Rate Dynamics during Meditation: An ECG Based Study with Multi-Fractality and Visibility Graph. Frontiers in Physiology. 7, (2016). https://doi.org/https://doi.org/10.3389/fphys.2016.00044

5. Boccaletti, S., Latora, V., Moreno, Y., Chavez, M., Hwang, D.-U.: Complex networks: Structure and dynamics. Physics Reports. 424, 175-308 (2006). https://doi.org/https://doi.org/10.1016/j.physrep.2005.10.009

6. Brú, A., Gómez-Castro, D., Nuño, J.C.: Visibility to discern local from nonlocal dynamic processes. Physica A: Statistical Mechanics and its Applications. 471, 718-723 (2017). https://doi.org/https://doi.org/10.1016/j.physa.2016.12.078

7. Chen, G., Wang, X., Li, X.: Fundamentals of Complex Networks: Models, Structures and Dynamics. Willey (2015)

8. Ding, D., Tang, Z., Wang, Y., Ji, Z.: Synchronization of nonlinearly coupled complex networks: Distributed impulsive method. Chaos, Solitons & Fractals. 133, (2020). https://doi.org/https://doi.org/10.1016/j.chaos.2020.109620

9. Dörfler, F., Francesco, B.: Synchronization in complex networks of phase oscillators: A survey. Automatica. 50, Pages 1539-1564 (2014). https://doi.org/https://doi.org/10.1016/j.automatica.2014.04.012

10. Estrada, E., Knight, P.A.: A First Course in Network Theory. Oxford University Press (2025)

11. Gonçalves, B.A., Carpi, L., Rosso, O.A., Ravetti, M.G.: Time series characterization via horizontal visibility graph and Information Theory. Physica A: Statistical Mechanics and its Applications. 464, 93-102 (2016). https://doi.org/https://doi.org/10.1016/j.physa.2016.07.063

12. Halvin, S., Kenett, D.Y., Ben-Jacob, E., Bunde, A., Cohen, R., Hermann, H., Kantelhardt, J.W., Kertész, J., Kirkpatrick, S., Kurths, J., Portugali, J., Solomon, S.: Challenges in network science: Applications to infrastructures, climate, social systems and economics. The European Physical Journal Special Topics. 273-293 (2012). https://doi.org/https://doi.org/10.1140/epjst/e2012-01695-x

13. Hou, F.Z., Li, F.W., Wang, J., Yan, F.R.: Visibility graph analysis of very short-term heart rate variability during sleep. Physica A: Statistical Mechanics and its Applications. 458, 140-145 (2016). https://doi.org/https://doi.org/10.1016/j.physa.2016.03.086

14. Jiruska, P., de Curtis, M., Jefferys, J.G.R., Schevon, C.A., Schiff, S.J., Schindler, K.: Synchronization and desynchronization in epilepsy: controversies and hypotheses. The Journal of Physiology. 591, 787-797 (2013). https://doi.org/https://doi.org/10.1113/jphysiol.2012.239590

15. Lacasa, L., Luque, B., Ballesteros, F., Luque, J., Nuño, J.C.: From time series to complex networks: The visibility graph. Applied Mathematics. 105, 4972-4975 (2008). https://doi.org/https://doi.org/10.1073/pnas.0709247105

16. Lacasa, L., Luque, B., Luque, J., Nuño, J.C.: The visibility graph: A new method for estimating the Hurst exponent of fractional Brownian motion. Europhysics Letters. 86, (2009). https://doi.org/10.1209/0295-5075/86/30001
17. Lacasa, L., Nicosia, V., Latora, V.: Network structure of multivariate time series. Scientific Reports. 5, (2015). https://doi.org/https://doi.org/10.1038/srep15508

18. Lacasa, L., Nuñez, A., Roldán, É., Parrondo, J.M.R., Luque, B.: Time series irreversibility: a visibility graph approach. The European Physical Journal B. 85, (2012). https://doi.org/https://doi.org/10.1140/epjb/e2012-20809-8

19. Li, R., Mayseless, N., Balters, S., Reiss, A.L.: Dynamic inter-brain synchrony in real-life inter-personal cooperation: A functional near-infrared spectroscopy hyperscanning study. NeuroImage. 238, (2021). https://doi.org/https://doi.org/10.1016/j.neuroimage.2021.118263

20. Long, X., Fonseca, P., Aarts, R.M., Haakma, R., Foussier, J.: Modeling cardiorespiratory interaction during human sleep with complex networks. Applied Physics Letters. 105, (2014). https://doi.org/https://doi.org/10.1063/1.4902026

21. Luque, B., Ballesteros, F., Núñez, A.M., Robledo, A.: Quasiperiodic Graphs: Structural Design, Scaling and Entropic Properties. Journal of Nonlinear Science. 23, 335-342 (2013). https://doi.org/https://doi.org/10.1007/s00332-012-9153-2

22. Luque, B., Lacasa, L., Ballesteros, F., Robledo, A.: Feigenbaum Graphs: A Complex Network Perspective of Chaos. PLOS ONE. 6, 1-8 (2011). https://doi.org/https://doi.org/10.1371/journal.pone.0022411

23. Mira- Iglesias, A., Conejero, J.A., Navarro-Pardo, E.: Natural visibility graphs for diagnosing attention deficit hyperactivity disorder (ADHD). Electronic Notes in Discrete Mathematics. 54, 337-342 (2016). https://doi.org/https://doi.org/10.1016/j.endm.2016.09.058

24. Newman, M.: Networks. Oxford University Press (2018)

25. Osipov, G.V., Kurths, J., Zhou, C.: Synchronization in Oscillatory Networks. (2007)

26. Pecora, L.M., Carroll, T.L.: Master Stability Functions for Synchronized Coupled Systems. Physical Review Letters. 80, (1998)

27. Pikovsky, A., Rosenblum, M., Kurths, J.: Synchronization. Cambridge University Press (2001)

28. Sannino, S., Stramaglia, S., Lacasa, L., Marinazzo, D.: Visibility graphs for fMRI data: Multiplex temporal graphs and their modulations across resting-state networks. Network Neuroscience. 1, 208-221 (2017). https://doi.org/https://doi.org/10.1162/NETN_a_00012

29. Tang, Y., Qian, F., Gao, H., Kurths, J.: Synchronization in complex networks and its application – A survey of recent advances and challenges. Annual Reviews in Control. 38, 184-198 (2014). https://doi.org/https://doi.org/10.1016/j.arcontrol.2014.09.003

30. Wang, X., Zhang, Y., He, Y., Lu, K., Hao, N.: Dynamic Inter-Brain Networks Correspond with Specific Communication Behaviors: Using Functional Near-Infrared Spectroscopy Hyperscanning During Creative and Non-creative Communication. Frontiers in Human Neuroscience. 16, (2022). https://doi.org/https://doi.org/10.3389/fnhum.2022.907332