## Abstract

We investigated whether functional brain networks are abnormally organized in Alzheimer's disease (AD). To this end, graph theoretical analysis was applied to matrices of functional connectivity of beta band–filtered electroencephalography (EEG) channels, in 15 Alzheimer patients and 13 control subjects. Correlations between all pairwise combinations of EEG channels were determined with the synchronization likelihood. The resulting synchronization matrices were converted to graphs by applying a threshold, and cluster coefficients and path lengths were computed as a function of threshold or as a function of degree *K*. For a wide range of thresholds, the characteristic path length *L* was significantly longer in the Alzheimer patients, whereas the cluster coefficient *C* showed no significant changes. This pattern was still present when *L* and *C* were computed as a function of *K*. A longer path length with a relatively preserved cluster coefficient suggests a loss of complexity and a less optimal organization. The present study provides further support for the presence of “small-world” features in functional brain networks and demonstrates that AD is characterized by a loss of small-world network characteristics. Graph theoretical analysis may be a useful approach to study the complexity of patterns of interrelations between EEG channels.

## Introduction

According to Delbeuck and others (2003), cognitive dysfunction in Alzheimer's disease (AD) could be due, at least in part, to a functional disconnection between distant brain areas. One possibility to examine this hypothesis is to study the correlations between signals of brain activity (EEG, magnetoencephalography [MEG] functional magnetic resonance imaging [fMRI] blood oxygen level dependent [BOLD]) recorded from different areas. The underlying assumption is that such correlations reflect, at least in part, functional interactions between different brain areas. The concept of statistical interdependencies between signals of brain activity as a tentative index of functional interactions is referred to as “functional connectivity” (for a review, see Lee and others 2003).

The synchronization likelihood (SL) is a recently introduced measure of statistical interdependencies within a dynamical systems framework (Stam and Van Dijk 2002). With this measure, a loss of upper alpha-, beta-, and gamma-band synchronization could be demonstrated in AD, both during a no-task state as well as during a working memory task (Stam and others 2002, 2003; Babiloni and others 2004; Pijnenburg and others 2004). In these studies, the 13- to 30-Hz beta band showed the most consistent abnormalities.

Although there seems to be growing consensus with respect to the loss of functional connectivity in AD, it remains unclear whether a decrease in the mean level of coupling is also associated with a change in the global organization of functional networks. Tononi and others (1998) have pointed out that optimal brain functioning requires a suitable balance between local specialization and global integration of brain activity. They indicated this optimal state as “complex” and proposed a neural complexity measure C_{N} that is sensitive to the optimal balance between segregation and integration (Tononi and others 1994). However, application of the neural complexity measure to fMRI, EEG, and MEG has not yet produced consistent results (Van Putten and Stam 2001; Burgess and others 2003; Van Cappellen van Walsum and others 2003; Branston and others 2005).

An alternative approach to the characterization of complex networks is the use of graph theory (Strogatz 2001; Sporns and others 2004). A graph is a basic representation of a network, which is essentially reduced to nodes (vertices) and connections (edges) (Fig. 1).

Graphs are characterized by a cluster coefficient *C* and a characteristic path length *L*, among other measures. The cluster coefficient is a measure of the local interconnectedness of the graph, whereas the path length is an indicator of its overall connectedness. Watts and Strogatz (1998) have shown that graphs with many local connections and a few random long distance connections are characterized by a high cluster coefficient and a short path length; such near-optimal networks are designated as “small-world” networks. Since then, many types of real networks have been shown to have small-world features (Strogatz 2001). Patterns of anatomical connectivity in neuronal networks are particularly characterized by high clustering and a small path length (Watts and Strogatz 1998). It has been suggested that a small world–like network architecture may be optimal for synchronizing neural activity between different brain regions (Lago-Fernandez and others 2000; Latora and Marchiori 2001; Barahona and Pecora 2002; Masuda and Aihara 2004). Networks of functional connectivity based upon recordings in animals, fMRI BOLD signals, or MEG recordings have also been shown to have small-world characteristics (Stephan and others 2000; Dodel and others 2002; Stam 2004; Eguiluz and others 2005; Salvador, Suckling, Coleman, and others 2005; Salvador, Suckling, Schwarzbauer, and Bullmore 2005).

In the present study, we intend to address the question whether functional brain networks in AD are characterized by a loss of small-world features such as a high cluster coefficient and a short path length.

## Materials and Methods

### Subjects

The study involved consecutive subjects referred to the Alzheimer Center at the VU University Medical Center. All subjects were studied according to a clinical protocol which involved history taking, physical and neurological examination, blood tests, MMSE, neuropsychological examination, magnetic resonance imaging of the brain, and a quantitative EEG. The final diagnosis was based upon a consensus meeting where all the available clinical data and the results of the ancillary investigations were considered. A diagnosis of probable AD was based upon the McKhann criteria (McKhann and others 1984). To overcome the lack of pathological verification, clinical diagnosis was monitored at regular intervals, and for studies, the unchanged diagnosis at 1 year was used. The same procedure has been used in other studies (De Leeuw and others 2004; Schoonenboom and others 2004).

The present study concerned 28 subjects, 15 with a diagnosis of probable AD (4 males, mean age 69.6 years, standard deviation [SD] 7.9, range 54–77) and 13 control subjects with only subjective memory complaints (subjective complaints [SC], 6 males, mean age 70.6 years, SD 7.7, range 57–78). Mean Mini Mental State Examination score of the Alzheimer patient group was 21.4 (SD 4.0, range 15–28); mean MMSE score of the SC subject group was 28.4 (SD 1.1, range 27–30).

### EEG Recording

EEGs were recorded in all subjects as part of the examination protocol. EEGs were recorded (against an average reference electrode) with an OSG digital EEG apparatus (Brainlab (R)) at the following positions of the 10–20 systems: Fp2, Fp1, F8, F7, F4, F3, A2, A1, T4, T3, C4, C3, T6, T5, P4, P3, O2, O1, Fz, Cz, and Pz. ECG was recorded in a separate channel. Electrode impedance was below 5 kΩ. Initial filter settings were: time constant = 1 s and low pass filter = 70 Hz. Sample frequency was 500 Hz and analog-digital precision 16 bit. EEGs were recorded in a sound attenuated, dimly lit room while patients sat in a slightly reclined chair. Care was taken by the EEG technicians to keep the patients awake during the whole recording. For the present analysis, 30 s of artifact-free data (containing no eye blinks, slow eye movements, excess muscle activity, electrocardiogram artifacts, etc.) were selected off-line. The EEG was down sampled to 125 Hz, resulting in time series of 4096 samples for further analysis. Digital, zero-phase shift filtering of the EEG in the beta band (13–30 Hz), computation of the SL, and the 2 graph theoretical measures (cluster coefficient and characteristic path length) were done off-line with the DIGEEGXP software written by one of the authors (C.S.). Graph theoretical analysis was based on the full 21 × 21 matrix of all possible pairwise combinations of electrodes.

### Computation of the SL

Correlations between all pairwise combinations of EEG channels were computed with the SL (Stam and Van Dijk 2002). Mathematical details can be found in the appendix to this paper; here we give a brief description. The SL is a general measure of the correlation or synchronization between 2 time series, which is sensitive to linear as well as nonlinear interdependencies. The SL ranges between P_{ref} (a small number close to 0) in the case of independent time series and 1 in the case of maximally synchronous signals. P_{ref} is a parameter that has to be set; in the present study, P_{ref} was set at 0.01. The basic principle of the SL is to divide each time series into a series of “patterns” (roughly, brief pieces of time series containing a few cycles of the dominant frequency) and to search for a recurrence of these patterns. The SL is then the chance that pattern recurrence in time series *X* coincides with pattern recurrence in time series *Y*; P_{ref} is the small but nonzero likelihood of coincident pattern recurrence in the case of independent time series. The end result of computing the SL for all pairwise combinations of channels is a square *N* × *N* matrix of size 21 (the number of EEG channels), where each entry *N _{i}*

_{,j}contains the value of the SL for the channels

*i*and

*j*.

### Computation of the Cluster Coefficient C and Characteristic Path Length L

The 1st step in applying graph theoretical analysis to synchronization matrices is to convert the *N* × *N* synchronization matrix into a binary graph. A binary graph is a network that consists of elements (also called “vertices”) and undirected connections between elements (called “edges”) (Fig. 1). Edges between vertices either exist or do not exist; they do not have graded values. The synchronization matrix can be converted to a graph by considering a threshold *T*. Because there is no unique way to choose *T*, we explored a whole range of values of *T*, 0.01 < *T* < 0.05, with increments of 0.001 and repeated the full analysis for each value of *T*. If the SL between a pair of channels *i* and *j* exceeds *T*, an edge is said to exist between *i* and *j*; otherwise no edge exists between *i* and *j*.

Once the synchronization matrix has been converted to a graph, the next step is to characterize the graph in terms of its cluster coefficient *C* and its characteristic path length *L*. A schematic explanation of graphs, cluster coefficients, and path lengths is given in Figure 1.

To compute the cluster coefficient of a certain vertex, we first determine to which other vertices it is directly connected; these other vertices (1 edge away) are called “neighbors.” Now the cluster coefficient is the ratio of all existing edges between the neighbors and the maximum possible number of edges between the neighbors; it ranges between 0 and 1. This cluster coefficient is computed for all vertices of the graph and then averaged. It is a measure for the tendency of network elements to form local clusters. The characteristic path length *L* is the average shortest path connecting any 2 vertices of the graph; the length of a path is indicated by the number of edges it contains. The path length *L* is an emergent property of the graph, which indicates how well its elements are integrated/interconnected.

When *C* and *L* are computed as a function of threshold *T*, the results might be influenced by differences in the mean level of synchronization between the 2 groups. Because the SL is expected to be significantly lower for Alzheimer patients than controls, for a given value of *T*, AD graphs will have fewer edges than controls graphs, and this will influence the differences in *C* and *L* between the 2 groups. To control for this effect, we repeated the analysis computing *C* and *L* as a function of degree *K*, which is the average number of edges per vertex. In this way, graphs in both groups are guaranteed to have the same number of edges so that any remaining differences in *C* and *L* between the groups reflect differences in graph organization.

The values of *C* and *L* as a function of degree *K* were compared with the theoretical values of *C* and *L* for ordered (*C* = 3/4, *L* = *N*/2*K*) and random (*C* = *K*/*N, L* = ln(*N*)/ln(*K*)) graphs. However, statistical comparisons should generally be between networks that have equal (or at least similar) degree sequences, as these are known to affect all kinds of network measures. Because the theoretical networks have Gaussian degree distributions and may thus not provide valid controls for the experimental networks in the present study, which may have some other degree distribution, we also generated random and ordered control networks following the procedure described by Sporns and Zwi (2004) and Milo and others (2002) which preserve the degree distribution exactly. For a *K* value of 3, for each EEG 20 random and 20 ordered networks were generated, and the mean *C* and *L* were calculated.

### Statistical Analysis

Statistical analysis consisted of independent samples *t*-tests and linear regression of the plots of *C* and *L* as a function of threshold. In order to investigate correlations between changes in topological parameters with cognitive measures, we calculated Pearson's correlation coefficient between MMSE scores (as a measure of cognitive function) and both cluster coefficient and path lengths.

## Results

As can be seen in Figure 2, the synchronization matrices of both groups show a complex but nonetheless rather similar pattern, with various regions of high (darker) and low (lighter) levels of synchronization. For instance, the dark region in the upper left corner corresponds to high levels of synchronization between prefrontal, frontal, and frontolateral channels.

Overall, the beta-band synchronization was lower in the Alzheimer group (main effect of group, *F*_{1,26} = 4.656, *P* = 0.040).

Figure 3 shows the graphs corresponding to the mean synchronization matrices of Figure 2 using a threshold *T* = 0.029. The graphs for both groups show a similar complex network, consisting of frontal and parieto-temporo-occipital components, linked by long distance connections (Fig. 3B,C).

Compared with the AD group, the graph of the control group has a larger number of edges between the central, temporal, and frontal regions (Fig. 3D). The graphs shown in Figure 3 represent group averages and serve primarily to illustrate the main patterns. For the actual analysis, the conversion of the synchronization matrix to a graph was done for each subject separately, and the averaging was done over the individual values of *C* and *L* as a function of the threshold *T*.

The mean cluster coefficient *C* as a function of threshold for the 2 groups is shown in Figure 4A.

For a low value of the threshold, the corresponding graphs are almost fully connected with edges between almost all vertices yielding a corresponding *C* close to 1 (for *T* = 0, *C* is expected to be 1). For increasing values of the threshold, more and more edges will be lost (providing the corresponding value of SL < *T*), and the cluster coefficient starts to decrease. Over the whole range of threshold values investigated (0.010–0.050), *C* of the control group is slightly higher than *C* of the Alzheimer group. However, due to the large variance, which increases with higher values of *T*, there are no consistent statistical differences between the 2 groups.

In contrast, the characteristic path length *L* does show clear differences between the groups (Fig. 4B). In the curves of Figure 4B, several patterns can be discerned. For small values of *T* (0.010–0.019), the path length increases almost linearly with the threshold. For increasing values of *T*, more and more edges will drop out, increasing the average path length between randomly chosen vertices. For an intermediate range of values of *T* (0.020–0.032), the path length is significantly larger for the Alzheimer group compared with the control group; the most significant difference is found for *T* = 0.029 (*t*-test, *P* = 0.0049). For further increases of *T*, the path length starts to level off; this phenomenon occurs earlier in the Alzheimer group than in the control group. This can be explained owing to the fact that for high values of *T* some of the vertices (EEG channels) become disconnected from the graph (splitting off); the resulting graph will be smaller, which will limit the further growth of *L*. For very high values of *T* (*T* > 0.043), the path length actually decreases with increasing *T*. This is due to the fact that the graph will be divided into 2 or more sub graphs (fragmentation). Because these sub graphs will be much smaller than the original full graph, the corresponding mean *L* will decrease. This graph fragmentation occurs earlier in the AD group than the control group, resulting in significantly larger path length of the controls for very high values of *T*.

Results of the analysis of *C* as a function of *K* are shown in Figure 5A.

As expected, *C* increases with *K*. No significant differences between the 2 groups are present. Comparisons with the theoretical values of *C* for ordered (*C* = 3/4) and random graphs with size *N* (*C* = *K*/*N*) (Fig. 5A) and with the constructed random and ordered matrices that preserve the degree sequences of their experimental counterparts (Fig. 5C) show that *C* for EEG data is intermediate between ordered and random graphs. Results of the analysis of *L* as a function of *K* are shown in Figure 5B. Here *L* decreases as *K* increases because more edges allow shorter possible paths. *L* is significantly longer in AD patients for 2.85 < *K* < 3.15. Comparison with theoretical values of *L* for ordered and random graphs (ordered graph: *L* = *N*/2*K*, random graph: *L* = ln(*N*)/ln(*K*)) shows that the path length of the EEG data is very short and smaller than that of a random graph for *K* < 3.4 (Fig. 5B). However, for the constructed random and ordered matrices, *L* of the EEG is lower than *L* of ordered networks and close to (but not smaller than!) *L* of random networks; also *L* of the AD group is significantly longer than *L* of the control group (Fig. 5D).

Pearson correlation coefficient between MMSE scores and path length was significant for the combined AD and control subject group: *r* = −5.91 (*P* = 0.01), but not for the AD group alone: *r* = −4.05 (*P* = 0.13). Correlations between MMSE scores and cluster coefficient were not significant for the combined AD and control subjects: *r* = 0.28 (*P* = 0.15) or the AD group: *r* = 0.22 (*P* = 0.42) (Fig. 6).

## Discussion

The principal finding of the present study is that changes in EEG beta-band functional connectivity display a loss of small-world network characteristics. We showed that AD was characterized by a longer characteristic path length with relative sparing of the local clustering. Functional connectivity matrices of beta-band synchronization were converted to graphs and analyzed in terms of cluster coefficients *C* and characteristic path length *L* for a range of thresholds. This approach is quite general and could also be applied to connectivity matrices based upon reconstructed sources in future studies. The main purpose of this analysis was to characterize the whole network in terms of local and global integration and to determine which aspect might be affected most in AD.

For a whole range of threshold values, the cluster coefficient showed a nonsignificant trend to lower values in AD patients (Fig. 4A), implying that the local connectedness of networks in AD is relatively spared. The sparing of the cluster coefficient in AD means that if a channel A is strongly correlated to 2 other channels B and *C*, then the likelihood that B and *C* will also be strongly coupled is not different from controls. An interesting finding was the different dependence of characteristic path length on threshold in Alzheimer patients and controls (Fig. 4B). Here the range of threshold values higher than 0.035 should be considered with caution because in this range, the size of the graph is no longer fixed (due to splitting off) and the graph may even split into sub graphs. For threshold values lower than 0.035, these problems do not exist, and a proper interpretation of *L* as the average shortest distance between any 2 vertices is valid. In this range, *L* was significantly longer for AD patients than for controls. Short path lengths have been shown to promote effective interactions between and across cortical regions (Sporns and Zwi 2004). Interactions between interconnected areas of the brain are believed to form the basis of cognitive processes (Pastor and others 2000; Friston 2002; Horwitz 2003). The significant correlation between path length and cognitive functioning (as measured with MMSE) in the present study is consistent with this notion and provides further evidence for the concept of AD as a disconnection syndrome. Figure 3D shows that the graph of the SL matrices in the AD group contains a lower number of edges between temporal, frontal, and specifically posterior central regions compared with the control group. Although the EEG has a poor spatial resolution and possible topographical implications of these findings must be interpreted with extreme caution, these areas are largely consistent with previously detected functional connectivity maps of the conscious resting state involving (among others) the temporal lobe, the posterior cingulate gyrus, and medial and lateral prefrontal cortex (Greicius and others 2003, 2004).

The finding of a lower level of beta-band synchronization is in agreement with several previous studies using the SL (Stam and others 2002, 2003; Babiloni and others 2004; Pijnenburg and others 2004). Claus and others (1998) also stressed the importance of the beta band in AD in the context of prognosis. The fact that the beta band shows fairly consistent changes in relatively mild AD was our reason for concentrating on this frequency band in the present study. The method described in the present study can now be used to analyze small-world characteristics in other bands, both during active and resting states.

A possible criticism of the graph theoretical results might be that they simply reflect the lower level of synchronization in the AD group. For a given value of *T*, graphs of AD patients are expected to have fewer edges, and this will result in a lower *C* and *L*. To determine whether the longer path length in the AD group reflects a true difference in organization and not simply a lower mean level of synchronization, we repeated the analysis by computing *C* and *L* as a function of *K*. In this way, graphs in both groups had equal numbers of edges, and any influence of differences in mean SL was eliminated. This analysis showed that AD patients still had a significantly longer *L* compared with controls, with no significant difference in *C*. Thus, the longer path length in AD patients cannot be ascribed to differences in mean level of SL and reflects a true abnormality in the organization of functional networks in this disorder.

A comparison of *C* and *L* of the EEG data with theoretical values of *C* and *L* for ordered networks (formula's taken from Watts and Strogatz 1998) showed that *L* of the EEG graphs was very small and even smaller than that of random graphs (Fig. 5B), whereas *C* of the EEG graphs was intermediate between ordered and random graphs (Fig. 5A). For a small-world graph, *L* should be close to that of a random graph, whereas *C* should be much higher than that of a random graph, and both would be expected to lie between the values of these 2 graphs. The reason the experimental *L* is smaller than that of random graphs may be because the experimental data do not have a Gaussian degree distribution. Therefore, we also generated random and ordered control networks following the procedure described by Sporns and Zwi (2004) and Milo and others (2002), which preserve the degree distribution exactly. *L* for the experimental data was close to (but not smaller than!) that of the random graph, whereas *C* was higher than that of the random graph, demonstrating that the pattern is still consistent with a small-world configuration. These results are comparable with those obtained for high- and low-frequency bands in a previous MEG study (Stam 2004).

In neuroscience, graph theoretical analysis has mainly been applied to the study of anatomical networks (Hilgetag and others 2000; Strogatz 2001; Sporns and others 2004). Different types of networks have been shown to be characterized by a relatively high cluster coefficient and a short path length, corresponding to the notion of small-world networks (Watts and Strogatz 1998; Sporns and others 2004). When networks are evolved while selecting for the highest complexity (defined as an optimal balance between local specialization and global integration), the resulting networks typically have the characteristic “small-world properties”; this suggests that small-world features and a high neuronal complexity C_{N} are in fact closely associated (Sporns and others 2000). There are indications that small-world networks and the related scale-free networks represent an optimal organization in terms of low “wiring costs,” local independence, and global integration. Modeling studies have shown that neural networks with a small-world configuration facilitate synchronization between distant neurons and efficient information processing (Lago-Fernandez and others 2000; Masuda and Aihara 2004). Stephan and others (2000) have shown that graph analysis can be applied equally well to patterns of functional and anatomical connectivity; in both cases a typical small-world network was revealed. In agreement with this, analysis of correlation matrices determined from fMRI BOLD signals has shown typical small-world patterns (Dodel and others 2002; Eguiluz and others 2005; Salvador, Suckling, Coleman, and others 2005; Salvador, Suckling, Schwarzbauer, and Bullmore 2005). Finally, in a study dealing with MEG recordings from healthy subjects, graph analysis of synchronization matrices revealed small-world patterns in low- and high-frequency bands (Stam 2004). When *C* and *L* are expressed as ratios of *C* and *L* of random graphs, the results of the present study are quite comparable with those of several previous studies (Table 1).

C/C_{random} | L/L_{random} | ||

Present study | AD | 1.6 | 1.12* |

Controls | 1.58 | 1.07 | |

Stam (2004) | Healthy subjects | 1.89 | 1.19 |

Salvador, Suckling, Coleman and others (2005) | Healthy subjects | 2.08 | 1.09 |

Hilgetag and others (2000) | Macaque visual cortex | 1.85 | 1.02 |

Cat whole cortex | 1.99 | 1.07 |

C/C_{random} | L/L_{random} | ||

Present study | AD | 1.6 | 1.12* |

Controls | 1.58 | 1.07 | |

Stam (2004) | Healthy subjects | 1.89 | 1.19 |

Salvador, Suckling, Coleman and others (2005) | Healthy subjects | 2.08 | 1.09 |

Hilgetag and others (2000) | Macaque visual cortex | 1.85 | 1.02 |

Cat whole cortex | 1.99 | 1.07 |

For the present study, data are shown for *K* = 3.0 for the generated control networks preserving *N, K*, and degree sequences; for the Stam (2004) study, data are shown for the gamma band and *K* = 20.

Significant difference between AD patients and controls, *P* = 0.046 (2-tailed *t*-test).

The present study provides further support for the presence of small-world features in functional networks in the brain. Furthermore, this study shows for the 1st time that pathological networks in AD may be less small world–like than normal brain networks. Alzheimer patients have significantly longer path lengths of their EEG graphs, even after correcting for differences in the mean level of synchronization, which suggests a disruption in effective interactions between and across cortical regions and provides further support for the concept of AD as a disconnection syndrome. Graph theoretical analysis can reveal abnormal patterns of organization of functional connectivity. This approach may be useful not only in degenerative dementia's but also in other disorders such as schizophrenia where abnormal functional connectivity plays a role (Friston 1999; Breakspear and others 2003).

### Appendix: Mathematical Details of Computation of SL

This appendix is based upon Posthuma and others (2005). The SL is a measure of the “generalized synchronization” between 2 dynamical systems *X* and *Y* (Stam and Van Dijk 2002). Generalized synchronization (Rulkov and others 1995) that exists between *X* and *Y* of the state of the response system is a function of the driver system: *Y* = F(*X*). The 1st step in the computation of the SL is to convert the time series *x _{i}* and

*y*recorded from

_{i}*X*and

*Y*as a series of state space vectors using the method of time delay embedding (Takens 1981):

*L*is the time lag and

*m*the embedding dimension. From a time series of

*N*samples,

*N*− (

*m*×

*L*) vectors can be reconstructed. State space vectors

*Y*are reconstructed in the same way.

_{i}SL is defined as the conditional likelihood that the distance between *Y*_{i} and *Y*_{j} will be smaller than a cutoff distance *r*_{y}, given that the distance between *X*_{i} and *X*_{j} is smaller than a cutoff distance *r*_{x}. In the case of maximal synchronization, this likelihood is 1; in the case of independent systems, it is a small, but nonzero number, namely, P_{ref}. This small number is the likelihood that 2 randomly chosen vectors *Y* (or *X*) will be closer than the cutoff distance *r*. In practice, the cutoff distance is chosen such that the likelihood of random vectors being close is fixed at P_{ref}, which is chosen the same for *X* and *Y*. To understand how P_{ref} is used to fix *r*_{x} and *r*_{y}, we first consider the correlation integral:

_{r}is the likelihood that 2 randomly chosen vectors

*X*will be closer than

*r*. The vertical bars represent the Euclidean distance between the vectors.

*N*is the number of vectors,

*w*is the Theiler correction for autocorrelation (Theiler 1986), and θ is the Heaviside function: θ(

*X*) = 0 if

*X*≥ 0 and θ(

*X*) = 1 if

*X*< 0. Now,

*r*

_{x}is chosen such that

*Cr*

_{x}= P

_{ref}, and

*r*

_{y}is chosen such that

*Cr*

_{y}= P

_{ref}. The SL between

*X*and

*Y*can now be formally defined as:

*X*and

*Y*(SL

_{XY}= SL

_{YX}). In equation (3), the averaging is done over all

*i*and

*j*; by doing the averaging only over

*j*, SL can be computed as a function of time

*i*. From equation (3), it can be seen that in the case of complete synchronization SL =1; in the case of complete independence SL = P

_{ref}. In the case of intermediate levels of synchronization P

_{ref}< SL < 1.

In the present study, the following parameters were used: *L* = 10, *M* = 10, and P_{ref} = 0.01.