Introduction

Measures of brain activity are good candidates to serve as endophenotypes for behavior and behavioral disorders, since most—if not all—behavior is governed by the brain. For this reason, many studies investigated the usefulness of electrophysiological measures to elucidate the pathway from gene to behavior. For example, electrophysiological measures of brain function such as the P300 have been suggested as an endophenotype of alcoholism (e.g., Porjesz and Begleiter 1990) and have been used in the identification of chromosomal regions in alcoholism (Williams et al. 1999). Another example is Frontal Asymmetry of resting state EEG as an endophenotype of depression and anxiety disorders (Anokhin et al. 2006; Coan and Allen 2004; Coan 2003; Smit et al. 2007).

These EEG endophenotypes often reflect the response to a specific class of stimuli or the activity in a particular brain region. To capture brain activity underlying behavioral traits it may be necessary, however, to additionally focus on the capacity of the brain for dynamic interaction between different brain regions. It has been argued that some disease states are linked to dysfunctional connectivity rather than dysfunctions of particular brain areas. Autism, for example, may be linked to increased local frontal brain interconnectivity with increased frontal white matter volume (Herbert et al. 2004), but reduced long range connectivity between the frontal lobe and the remaining brain areas (Coben et al. 2008; Courchesne and Pierce 2005; Koshino et al. 2005). In the non-clinical field, Thatcher et al. (2005) have shown that measures of connectivity outperform other psychophysiological measures in explaining individual differences in IQ.

Due to their high temporal resolution, the electroencephalogram (EEG) and magnetoencephalogram (MEG) may be particularly useful for measuring functional connectivity of the brain. Generally, coherence has been the measure of choice to provide information on the overall connectedness between two brain regions whose activity is often represented by two separate EEG or MEG signals (Thatcher et al. 1987, 2008). Coherence captures the linear statistical dependencies in the frequency domain between pairs of signals. The brain, however, may perhaps better be viewed as a collection of dynamically interacting regions (Friston 2000; Pereda et al. 2005). To measure this type of connectivity, a measure has been proposed that captures not just linear, but also non-linear interactions (Stam and van Dijk 2002). Synchronization Likelihood measures the coupling strength between two systems (e.g., brain regions), each of which is represented by a signal. The SL, like coherence, captures overall functional connectivity but additionally detects non-linear coupling, does not show spurious connectivity between bandpass filtered signals, is insensitive to fluctuations in power, and robust to non-stationarity of the EEG (Stam and van Dijk 2002).

Recently, it has become clear that the pattern of functional brain connectivity, as measured by the SL has a very specific structure. This structure is revealed by an analysis strategy that is based on graph theory and was initially proposed by Watts and Strogatz (1998). In their landmark article, the authors reduced networks to their mathematical essence: nodes and their connections. All types of these simplified networks (called graphs) can be functionally described by two parameters. The first describes the amount of local connectivity: clustering coefficient CC. It takes a value between 0 and 1 indicating the proportion of neighboring nodes that are interconnected amongst each other. The second parameter describes global interconnectedness and is called the average path length L. It is a value simply indicating the average number of steps required to go from a node to all others taking the shortest route. It has been argued that these parameters (and their derivations) describe non-trivial aspects of connection patterns. Clustering coefficient CC is a measure for low building cost and robustness against perturbations. Path length L represents the ability of the system to integrate information across distant sources (Achard et al. 2006; Stam 2004).

The graph theoretical approach allows one to describe networks along the dimension of ordered to random. Watts and Strogatz (1998) showed that ordered graphs (lattices, regular networks) are characterized by high CC and long L. Random graphs, i.e., with all connections randomly reshuffled between the nodes, have short L but low CC. By starting from ordered graphs and slowly randomly reconnecting single connections with a very small rewiring probability p, Watts and Strogatz showed that the path length L drops surprisingly quickly maintaining a high clustering coefficient CC. These types of networks combine the best of both worlds (low building costs, high integration), and are called Small-World networks. Many networks have been proven to have a small-world topology, such as human society (Milgram 1967), power grids, the world wide web, and many types of collaboration networks (Newman 2003). An increasing number of articles have shown that the resting human brain, too, has a small world topology (see Bassett and Bullmore 2006). MEG studies using graph theory to describe connectivity patterns in different frequency bands have revealed a clear small world topology (Stam 2004; Bassett et al. 2006). A small world topology has also been found in the classical frequency bands of the EEG (Micheloyannis et al. 2006).

All three connectivity measures (SL, CC and L) show large individual differences that appear partly driven by genetic factors (Posthuma et al. 2005; Smit et al. 2008). This indicates that these measures could serve as endophenotypes for many aspects of normal and abnormal brain function. In this article, we will expand on previous studies in two ways. First, we will investigate whether the three connectivity parameters reflect a single source of genetic variation. That is, are the graph parameters influenced by overlapping sets of genes? Second, we will investigate the longitudinal stability and genetic stability spanning a period from adolescence into young adulthood. We invited twins and siblings into the EEG laboratory for resting EEG measurement as part of longitudinal studies into the genetics of brain and behavior by the Netherlands Twin Registry (NTR; Posthuma et al. 2001, 2005; van Beijsterveldt et al. 1996). The current study combined longitudinal EEG data in twins measured at age 16, and returning at age 18 and 25 with EEG data in two additional cohorts of young (ca. 25 years) and middle-aged (ca. 50 years) twins and siblings.

Methods

Subjects

The sample for this study was derived from an ongoing twin family study on cognition in twins and family members from the Netherlands Twin Registry (Posthuma et al. 2001; van Beijsterveldt et al. 1996). Twins and siblings were invited for detailed psychophysiological study in the laboratory including eyes closed resting EEG and several experimental EEG task designs. The adolescent sample consisted of 213 twin pairs (75 MZ and 138 DZ, including 48 opposite sex pairs) with average age of 16.05 years (SD = .55 years, range 14.84–17.99) (van Beijsterveldt et al. 1996). 180 of the twin pairs revisited the laboratory on average 1.52 years later (SD = .10 years, range 1.29–2.50), and 53 twin pairs on average 5.85 years after the last session (SD = .75). The latter twins were part of a larger adult EEG sample that consisted of 760 subjects (twins as well as siblings) from 309 families divided into two age cohorts based on the age of the twins: a younger cohort of 159 families (mean age of twins = 25.8 years, SD = 2.9, range 18.75–33.9) and a middle-aged cohort of 150 families (mean age of the twins = 49.4 years, SD = 6.8, range 36.0–71.0). This sample consisted of 141 families with an MZ twin, 168 families with a DZ twin (62 of which were opposite sex twins). Participating families consisted of one to seven siblings (including twins). On average, 2.50 participants per family participated (Smit et al. 2005, 2008). Informed consent was obtained in writing for the EEG study. The study received approval from the appropriate ethical committees.

EEG recording

The experimental protocol for background EEG registration has been described in detail elsewhere (van Beijsterveldt et al. 1996; Posthuma et al. 2001; Smit et al. 2005). A brief description is given here. Resting EEG for the adolescent sample was measured using a Nihon Koden 4414A1K electroencephalograph and Ag/AgCl electrodes in a lycra Electro-cap in positions Fp1, Fp2, F3, F4, F7, F8, C3, C4, P7, P3, P4, P8, O1, and O2. Tin reference electrodes were attached to the ears and linked using the high input impedance pre-amplifier method of Pivik et al. (1993). This method avoids artifacts due to imbalanced electrode impedance. The vertical electro-oculogram (EOG) was recorded bipolarly between two tin electrodes, affixed one cm below the right eye and one cm above the eyebrow of the right eye. The horizontal EOG was recorded with tin electrodes affixed 1 cm left from the left eye and 1 cm right from the right eye. Another tin electrode placed on the forehead served as a ground electrode.

Resting EEG for the adult sample was measured using either an amplifier developed by Twente Medical Systems (TMS; Enschede, The Netherlands) for 661 subjects (375 young, 286 middle-aged) or a NeuroScan SynAmps 5083 amplifier for 104 subjects (24 young, 80 middle-aged), both systems using Ag/AgCl electrodes mounted in a lycra Electro-cap including positions F1, F2, F3, F4, F7, F8, C3, C4, P3, P4, P7, P8, O1, and O2. Reference electrode was A1, but later digitally rereferenced to linked ears. Amplifier filter settings for TMS were a single order FIR bandpass filter with cutoff frequencies of 0.05 and 30.0 Hz. NeuroScan filter settings were a lowpass filter at 50.0 Hz. EOG was measured with Ag/AgCl electrodes on the same positions as the adolescent groups.

For all groups data were digitized at 250 Hz. Impedances of all electrodes were kept below 10 kΩ. Subjects were seated in a comfortable reclining chair in a dimly lit, sound attenuated, and electromagnetically shielded room. They were instructed to close their eyes, relax, but stay awake and minimize eye and body movement.

EEG data processing

Three subjects’ data were unavailable because of technical difficulties resulting in lost data. All available EEG was visually checked for bad channels such as absence of signal, hum, clipping, persistent muscle tone artifacts, and external noise. In addition, bad episodes were marked. Subjects without the full set of 14 leads were excluded. Next, the data were filtered using a Matlab FIR filter from 2 to 38 Hz (zero-phase shift, 6 dB roll off). To remove eye artifacts, we used the ICA filtering technique implemented in EEGLAB, an open source MATLAB toolbox (Delorme and Makeig 2004). Independent components were determined in the combined EEG + EOG data. EOG artifacts were then selected on the basis of showing a specific frontal loading (asymmetric loadings for horizontal eye movements and symmetric for vertical eye movements) and a high correlation with one of the EOG channels (r > .7). These components were then removed via backprojection of the component activations using only the remaining components.

Data were then cut into 8 s epochs of artifact-free signals. From these, 12 were randomly selected. Subjects with less than twelve epochs were excluded. In total, none of the 16-year-olds, 2 of the 18-year-olds, 116 of the adult cohort were excluded. Each epoch was filtered into alpha (7.0–13.0 Hz), and beta (15.0–25.0) activity.

Synchronization likelihood

If a signal s1 is in a certain state (to be defined below) at time i we may find a recurrence of that state at another time points j. Next, we look if the second signal s2 is in the same state at time points i and j, recording a hit if so. SL is defined as the proportion of hits to the total number of recurrences in the s1 and is thus a number between 0 and 1. Note that SL is found even when signals s1 and s2 are in different states at i and j, as long as s1 and s2 are self-similar at i and j.

More formally, the instantaneous state of an EEG signal was represented by m-dimensional state vectors embedded vector X i = {x i, x i + 1l, x i + 2l, …, x i + (m − 1)l} where l is the lag and m the embedding dimension. The elements of X i are m samples taken from the signal spaced l apart. The vector is taken to represent the state of the system at time i. Within the same signal recurrences are sought at times j that reflect a similar state: A threshold distance εx is chosen such that a fixed proportion (p ref) of comparisons are close enough to be considered in a similar state. Next, the same comparison is made for a different system Y at the same time points i and j and with the same value for p ref. Now the synchronization likelihood S i between X and Y at time i is defined as follows:

$$ S_{i} = {\frac{1}{N}}\sum\limits_{j} {\theta \left( {\varepsilon_{y} - \left| {Y_{i} - Y_{j} } \right|} \right)\theta \left( {\varepsilon_{x} - \left| {X_{i} - X_{j} } \right|} \right)} $$

where θ is the Heaviside step function returning 0 for all values <0 and 1 for values ≥0. N represents the number of recurrences in signal X, i.e.:

$$ \sum\limits_{j} {\theta \left( {\varepsilon_{y} - \left| {X_{i} - X_{j} } \right|} \right)} $$

Overall SL between X and Y is the average over all possible i. To withhold the system to compare X i and X j while they represent the same state, only values for j are considered that are at sufficient time distance. The value of this distance, W1, is the Theiler correction for autocorrelation (Theiler 1986). The value for |i − j| is upper bound to create a window (W1 < W2 < N) to sharpen the time resolution of S i . More details on SL calculation can be found in several other publications (Montez et al. 2006; Posthuma et al. 2005; Stam and van Dijk 2002). The parameter settings l, m, W1 and W2 were chosen based on the filter frequency settings (see Montez et al. 2006). W2 and pref had fixed values of W1 + 400 and 0.02. These values reflect similar choices from the previous literature (Ponten et al. 2007; Smit et al. 2008), but p ref was increased from 0.01 to 0.02 to reflect the lower number of recurrences possible in the shorter epoch length of 8 s. These parameter settings are not critical (e.g., see Smit et al. 2008; Stam and van Dijk 2002).

CC and L calculation

SL was computed between each pair of electrodes resulting in a (14, 14) matrix where the values on the diagonal will be ignored. A binary graph was formed by applying a threshold θ such that the average degree was fixed at K = 4, that is, each node had on average 4 connections. Clustering coefficient CC and average path length L were calculated following Watts and Strogatz (1998) and as explained in detail elsewhere (Smit et al. 2008). In addition, we followed Newman’s (2003) proposal to assign the value of +∞ to the path length involving unconnected nodes and use the harmonic mean to average the values obtained for path length L across the different nodes.

Previous results had shown that measurement reliability of SL in four 16 s epochs was relatively high (r > .95), but that of the graph parameters CC and L derived from these SL estimates was much lower (for CC: r < .62; for L: r < .82; Smit et al. 2008). We suspected that this discrepancy might be ameliorated by using more epochs of shorter period. Therefore, we used 12 epochs of 8 s recordings, which increases the number of graphs calculated from the EEG signals.

EEG power

All three connectivity measures SL, CC, and L are to some degree sensitive to signal-to-noise ratios. Since some subjects may be lacking certain oscillations that are ubiquitous in the population (e.g., the absence of alpha oscillations in low voltage subjects; Vogel 2000), they may show spuriously low scores on the connectivity parameters as well. We therefore measured power of the oscillations in the alpha and beta bands using the pmtm function in MATLAB for each channel and frequency band implementing the Thompson (1982) multitaper method for estimating power spectral density. The average across channels was used to represent average EEG power in a frequency band. EEG power was included in the multivariate analysis.

Genetic and environmental variance decomposition

For all statistical modeling the freely available software package Mx version 1.65a was used (Neale 2004).

Since DZ correlations were less than half the MZ correlations in nearly all variables and not significant otherwise we did not investigate the effects of common environment (see “Results” section). Instead, dominant genetic factors were used in all variance decomposition models. Decomposition of variance consisted of Cholesky decompositions into Genetic (Additive and Dominant genetic) and Unique Environmental variance of the observed scores of SL, CC, L, and power. A tetravariate Cholesky decomposition was used to investigate genetic and environmental covariation between the four measures: SL, CC, L, and EEG power (see Fig. 1). These were performed on each of the four age groups separately. Means were modeled as a function of age group and sex.

Fig. 1
figure 1

Path model describing the multivariate variance decomposition across variables for age group 16 repeated for the other three age groups. Variance in the observed variables SL, C, L, and Power was explained in a Cholesky decomposition of environmental variance by 4 factors (E1–E4) that were uncorrelated between individuals 1 and 2 (no arrows), and genetic variance by 4 additive genetic factors (A1–A4) that correlated 1.0 between MZ twin pairs and 0.5 between DZ twins and siblings, and by 4 dominant genetic factors (D1–D4) that correlated 1.0 between MZ twins and 0.25 between DZ twins and siblings. The model can be expanded with additional siblings

To examine longitudinal stability in subjects that participated at age 16, 18, and 25, a trivariate Cholesky decomposition was used for each phenotype separately (see Fig. 2). Data from age group 50 did not have any overlap in sample with the other three age groups, and were therefore not included in the longitudinal models. Means were modeled as a function of age cohort and sex.

Fig. 2
figure 2

Path model describing the trivariate model for longitudinal variance decomposition of SL. The model was repeated for C and L. Variance in the observed variable (i.e., SL) measured at three time points (16, 18, and 25) was explained in a Cholesky decomposition of environmental variance by factors E1 to E3 that were uncorrelated between individuals (no arrows), and genetic variation by 4 additive genetic factors (A1–A3) that correlated 1.0 between MZ twin pairs and 0.5 between DZ twins and siblings, and by 4 dominant genetic factors (D1–D3) that correlated 1.0 between MZ twins and 0.25 between DZ twins and siblings. The model can be expanded with additional siblings

To calculate phenotypic and genetic correlations (multivariate model) or phenotypic and genetic stability (longitudinal model) we used the estimates from the full ADE decomposition. Since genetic effects comprised both Additive genetic and Dominant genetic effects, we estimated an overall statistic that summarizes the broad sense genetic correlation (broad r G) by dividing the total genetic covariation between x and y by the total genetic variation:

$$ {\text{broad}}\,r_{G} = {\frac{{Cov_{A,xy} + Cov_{D,xy} }}{{\sqrt {V_{A,x} + V_{D,x} } \times \sqrt {V_{A,y} + V_{D,y} } }}} $$

where x and y represent two different variables or two different time points.

To test the effects of age on SL, CC, and L we fitted first (linear) and second order (quadratic) polynomial regressions by including age and age2 predictors in the means prediction in Structural Equation Models taking within-family correlational structure into account.

All testing was performed with an alpha level of 0.01.

Results

Monozygotic and dizygotic twin correlations were estimated in univariate saturated models (see Table 1). In all cases but two, DZ correlations were less than half the MZ correlations. In the two other cases, common environment was not significant. Therefore, dominant genetic factors were included in the multivariate and longitudinal variance decomposition models.

Table 1 Twin correlations from univariate models for connectivity parameters SL, CC, and L in age groups 16–50

Multivariate models at each age

Phenotypic and genetic correlations between the variables SL, CC, L, and EEG power were estimated in the full ADE models (Fig. 1) and shown in Fig. 3. The top panels show the phenotypic correlations between the connectivity measures SL, CC, and L for each of the age groups, bars indicate the 99% confidence intervals. SL and L showed moderate (beta band, .40 < r < .53) to strong (alpha band, .56 < r < .73) correlation, all highly significant. Genetic overlap was high in both the alpha (.73 < broad r G < .79 across age groups) and beta band (.55 < broad r G < .68 across age groups). In addition, both SL and L were moderately related to EEG power in the respective frequency bands in the alpha band (SL and power: .30 < r < .46; L and power: .36 < r < .52). The genetic correlation was of a similar magnitude (SL and EEG power: .33 < broad r G < .53; L and EEG power: .41 < broad r G < .64). In the beta band, the phenotypic and genetic correlations of SL and L with EEG power were either non-significant, less strong, and/or varied in direction depending on age.

Fig. 3
figure 3

Phenotypic and genetic correlations between connectivity parameters SL, CC, L, and EEG power from the multivariate statistical models. See text for the interpretation of the results

CC was a largely independent measure and did not show a consistent or strong relationship with both L and SL (alpha band: −.30 < r < .15; beta band: −.28 < r < .27). Genetic correlations of CC with SL and L were at best low (broad r G < .30) that did not reach significance in most cases. The single (anomalous) exception was a large and significant genetic correlation between CC and L at age 16 (broad r G = .84). In addition, there was no substantial relation to EEG power phenotypically (−.02 < r < .20) or genetically (.05 < broad r G < .34).

Longitudinal models and heritability

Figure 4 shows the changes of the graph and connectivity parameters across age. SL showed highly significant quadratic regression terms (alpha: χ2 = 138.7, df = 1, p < 10−31; beta: χ2 = 116.8, df = 1, p < 10−26). For CC, a linear fit was significant in the alpha band (χ2 = 17.45, df = 1, p < .0001), but no age regression was found for CC in the beta band. L showed highly significant quadratic polynomial regression coefficients in both the alpha (χ2 = 28.6, df = 1, p < 10−7) and beta band (χ2 = 18.2, df = 1, p < .0001). EEG power was also included to test for possible dependencies of the connectivity parameters on overall signal strength of the alpha and beta oscillations. EEG power showed no significant regression in the alpha band (χ2 = 6.2, df = 1, p = .011) but a significant 2nd order polynomial in the beta band (χ2 = 8.47, df = 1, p < .01). Given the differential developmental patterns in EEG power and the connectivity measures it is unlikely that power caused the age effects in connectivity.

Fig. 4
figure 4

Plots of the effect of age on the connectivity parameters. Left column for Alpha oscillations (7–13.0 Hz), right column for Beta oscillations (15–25 Hz). Plots show linear (1st order) and quadratic (2nd order) polynomial regressions when significant. Dashed lines indicate values for the graph parameters obtained by randomizing all connections in the graphs (CCrnd, L rnd) or by calculating connectivity between white-noise signals (SLrnd). L show an inverted U in both bands with a maximum at ~40 years. L was close but systematically higher than L rnd. CC showed only slight changes with age (alpha band) and was much higher than CCrnd at all ages and for all individuals. L ≈ L rnd and CC ≫ CCrnd indicates a small world connectivity pattern and can be found at all ages. Synchronization Likelihood (SL). Overall synchronicity in the signal develops similar to L in an inverted U. EEG power showed only developmental changes for beta

We calculated phenotypic and genetic stability from full ADE models. The upper panels of Fig. 5 show the stability of the connectivity parameters between 16 and 18, 16 and 25, and 18 and 25. In the alpha band, stability measures are generally high and significant. Stability was >.73 for SL, >.52 for CC, and >.59 for L. Beta band stability was lower ranging from about 0.0–0.6 and failed to reach significance for beta band CC between 16 and 25 and 18 and 25. There was some evidence for a decrease in stability with age as evidenced by the 99% confidence intervals. Significance of a difference may be inspected by determining whether the stability values fall out of each other’s CI. Phenotypic correlation for alpha band L from 16 to 25 was significantly lower than stability from 16 to 18, and marginally significantly lower than 18–25 (i.e., stability 16–25 was not included in the CI of stability 18–25, but stability 18–25 was included in the CI of stability 16–25). Beta band CC stability from 16 to 18 was significantly higher than stability across longer age ranges (18–25 and 16–25).

Fig. 5
figure 5

Longitudinal stability and genetic longitudinal stability between age groups 16, 18, and 25 for connectivity parameters SL, CC, and L. Scores in the alpha band showed all the hallmarks of a good endophenotype, i.e., moderate to high stability and high genetic stability. Beta band stability was much lower (especially clustering coefficient CC). Genetic correlations were not significantly different from unity, but in many cases also not significantly different from zero

The lower panels of Fig. 5 show the broad genetic stability (broad r G) across the age groups with longitudinal data available (viz., 16, 18, and 25). Genetic stability was high for alpha band parameters, significantly higher than zero and not significantly different from unity suggesting high genetic stability. But in some cases significance could not be established due to a lack of statistical power that accompanies genetic correlations between factors with low heritabilities; For example, even though the genetic stability of CC in the beta band between ages 16 and 25 was 0.88, it was not significantly different from zero because the underlying genetic factor for age 16 was only just significant (h 2 = 25%, see Table 2).

Table 2 Heritabilities and 99% confidence intervals for SL, CC, L, and EEG power estimated with the multivariate models

Discussion

We computed overall SL from matrices of SL values calculated between pairs of resting state EEG signals. In addition, we applied the graph theoretical approach to the SL matrices to derive the parameters CC and L. These three functional connectivity measures were then subjected to a series of genetic analyses that tested their viability as endophenotypes in genetic association and linkage research.

In the alpha band, all three parameters reflected stable and moderately to highly heritable traits. Alpha SL proved to be a highly stable measure across a 9-year period (.73 or higher) and highly heritable. In adolescence heritability estimates (82% at age 16, 77% at age 18) were about as high as the heritability of EEG power (van Beijsterveldt et al. 1996). Heritability decreased in adulthood (62% at age 25 and 53% at age 50) but genetic factors remained the main source of variance. CC and L were less stable across the 9-year period (between .52 and .73) and less heritable (37–63%). Importantly, however, genetic stability was high (>.81) from 16 to 25 years for all three connectivity parameters CC, L, and SL.

Beta band activity showed similar results, with genetic influences generally more modest than in the alpha band. Heritability of beta SL ranged between 40 and 75%, and heritability of L between 29 and 42% and CC between 25 and 40%. Temporal stability of the beta band connectivity parameters was lower than in the alpha band, and not significant for the CC parameter. As in the alpha band, genetic stability was moderate to high (mostly >.60) and not significantly different from unity. However, they showed wider confidence intervals and were often also not significantly different from zero. Lower heritability and genetic stability may in part reflect the lower overall power in the beta band—and therefore lower signal-to-noise ratios. We suspect that beta band connectivity measures may require much longer EEG recordings to sufficiently reduce measurement error.

Both alpha and beta band heritability of the graph parameters CC and L were similar to our previous findings from the adult subset of the current sample (Smit et al. 2008), although a straightforward comparison may not be entirely feasible due to differences in epoch selection (12 times 8 s epochs against 4 times 16 s epochs), and filtering (two wide against five narrow frequency bands). Also in keeping with previous results we found no evidence for shared environment on the connectivity measures, and DZ correlations were well below half the magnitude of the MZ correlation.

Correlations across the different connectivity parameters revealed that there is some interdependency between SL and L. The correlation between SL and L was similar in all age groups (ca. 0.65 for alpha and ca. 0.45 for beta), and largely caused by shared genetic factors (broad genetic correlation ca. 0.80 for alpha and 0.60 for beta). This could indicate that the average path length L is a parameter that is sensitive to the overall connectivity strength, even though graphs were constructed using variable thresholds, i.e., corrected for the overall level of SL. CC showed low correlations with both SL and L, and thus seems to reflect different aspects of the organization of the connectivity pattern. In addition, genetic correlations of CC with SL and L showed an overall pattern of being nonsignificant or low in magnitude suggesting that CC may index different genetic sources of variation than L and SL. Note that confidence intervals around these genetic correlations were fairly large, which may be an indication of a lack of statistical power instead of truly separate genetic sources of variation.

The connectivity parameters showed clear changes with age. Overall synchronization, for instance, showed an increase from adolescence to young adulthood, peaking at 40–50 years of age, and a decrease again into older age. Path length L showed a similar inverted U shape, a further indication of the interdependence of SL and L. It also indicates that the resting state functional network resembles a more randomly connected network in both adolescent and older age groups, and a more organized structure during young and middle-aged adulthood. Clustering coefficient CC did not show developmental changes, although a significant decrease with age was found for alpha band connectivity. Figure 4 also showed that L is quite near the L of a random graph (L rnd; dashed line in the figure) while CC was clearly higher than CCrnd. In addition, absolute values of CC were relatively stable over the years. Therefore, at all ages the brain showed the small-world network architecture characterized by L ≈ L rnd and CC≫CCrnd (Stam 2004). Recently, a study reported the development of SL, C, and L from childhood (ranging 8–12 years) to young adulthood (ranging 21–26 years; Micheloyannis et al. 2009). Some significant differences were reported between these two age groups, but these were restricted to the beta and gamma bands for resting state EEG. Interestingly, SL showed higher values in childhood, suggesting that the inverted-U shape of Fig. 4 may not continue into lower ages. CC also showed a decrease from childhood to adulthood, suggesting that this parameter is a marker of childhood development and stays relatively stable in adulthood.

It might be argued that the observed age effects have been caused by a lack of signal magnitude—and thus a change in measurement quality rather than a ‘true’ developmental change within the subjects. We investigated this issue by assuming that EEG power of alpha and beta oscillations may be an indication of the amount of signal available in the ongoing EEG. Developmental curves for EEG power were either absent or very unlike the clear polynomial regression coefficients found for SL and L. Therefore, the age related changes in parameters SL and L are unlikely to have been caused by a developmental change in EEG power. Likewise, it could be argued that individuals with low EEG power generally have graph parameters closer to the random situation than those with high power due to decreased signal-to-noise ratios. We inspected this by correlating EEG power with SL, CC, and L. The results indicated that in the alpha band power was moderately but clearly and positively related to both L (.36 < r < .52) and SL (.30 < r < .46) explaining on average 20 and 16% of the variation in L and SL, respectively. It was, however, largely unrelated to CC. In the beta band power correlated significantly only to L (r = .21) but failed to explain a substantial amount of variation (4%). Therefore, power may explain some of the individual variation connectivity parameters SL and L, but this effect is restricted in magnitude and to alpha band connectivity.

The current results were limited in several ways. First of all, we did not have a fully compatible set of electrodes between the age groups. For age groups 16 and 18, we had frontal leads Fp1 and Fp2 available. For age groups 25 and 50, F1 and F2 were available. This limits the interpretability of the results in some measures. However, the (genetic) stability between adolescent and adult age groups may serve as an indication that the substitution of leads did not alter the results much. Secondly, we have no certainty that volume conduction could have caused spurious synchronicity between adjacent electrode pairs, thus increasing overall SL and shortening path Length L. In addition, we have provided no indication of whether all age related changes did result from such a spurious effect. The development of a single, large neural generator in one particular area—for example, in the frontal brain—may lead similar effects of increased path length and decreased L through spurious conductivity effects. However, it should be noted that electrodes were sparse and many combinations of electrodes were close to or further apart than the 8 cm necessary to avoid a large influence of volume conduction (Nunez and Srinivasan 2006). Thirdly, many DZ correlations were well below half the magnitude of the MZ correlation. Although our sample was larger than any previous sample and the use of more epochs of shorter duration increased heritability compared to previous reports (Smit et al. 2008) confidence intervals around genetic correlations were wide. Larger samples and heritability estimates less influenced by measurement noise may perhaps ensure that the non-significant correlations between CC and other parameters hold.

A final limitation may lie in the uncertain nature of resting state condition. The resting state was chosen for its large history in the extant literature as well as its proven clinical status. Nevertheless, the uncontrolled free-flow of thought may yield interindividual differences in the brain network activation, thus confounding the current results: MZ twins may respond similarly to the task instructions by, for example, counting, or mentally creating grocery lists. However, recent studies on resting-state fMRI have examined the temporal fluctuation of the BOLD signal and extracted networks of coactivated voxels using Independent Components Analysis. The resulting resting-state networks are highly consistent between subjects—including the so-called default-mode network (Damoiseaux et al. 2006)—as well as highly consistent between task and resting conditions (Smith et al. 2009). These findings suggest that the resting state condition is not the source of interindividual differences and generalizes to the cognitively active brain. As a final note, it may be argued that interindividual differences in network activation cause behavioral differences in nature and content of thought.

Often cited criteria for endophenotypes are heritability, stability and theoretical meaningful ‘intermediateness’ of the endophenotype for behavioral outcomes (de Geus 2002; Gottesman and Gould 2003). Since many genetic association and linkage studies use subjects of different ages, an additional criterion is genetic stability, where the same set of genes is seen to influence the endophenotype across a wide age range. Taken together, our results suggest that the connectivity measures in the alpha band fulfill three of the criteria of an endophenotype: They are phenotypically stable and heritable traits, and the underlying genetic factors are stable. We further argue that connectivity measures are theoretically meaningful constructs to understand the brain phenotype. Small-world networks—which combine high clustering with low average path length—show resilience to perturbations and allow efficient information exchange (Achard et al. 2006). In addition, they combine modularity with global information integration (Bassett and Bullmore 2006). Moreover, a small-world topology may be a natural end-state of neuronal networks to evolve into (Siri et al. 2007).

In conclusion, we have shown that the connectivity parameters SL, C, and L all showed the hallmarks of a good endophenotype, at least in the alpha band, and may provide insight into how the brain develops—and thus how brain development may fail. Similar measures in the beta band may require prolonged EEG registrations and a larger number of epochs to obtain more reliable estimates of CC and L from the SL connectivity matrix.