Disruptions in cortico-subcortical covariance networks associated with anxiety in new-onset childhood epilepsy.

Anxiety disorders represent a prevalent psychiatric comorbidity in both adults and children with epilepsy for which the etiology remains controversial. Neurobiological contributions have been suggested, but only limited evidence suggests abnormal brain volumes particularly in children with epilepsy and anxiety. Since the brain develops in an organized fashion, covariance analyses between different brain regions can be investigated as a network and analyzed using graph theory methods. We examined 46 healthy children (HC) and youth with recent onset idiopathic epilepsies with (n = 24) and without (n = 62) anxiety disorders. Graph theory (GT) analyses based on the covariance between the volumes of 85 cortical/subcortical regions were investigated. Both groups with epilepsy demonstrated less inter-modular relationships in the synchronization of cortical/subcortical volumes compared to controls, with the epilepsy and anxiety group presenting the strongest modular organization. Frontal and occipital regions in non-anxious epilepsy, and areas throughout the brain in children with epilepsy and anxiety, showed the highest centrality compared to controls. Furthermore, most of the nodes correlating to amygdala volumes were subcortical structures, with the exception of the left insula and the right frontal pole, which presented high betweenness centrality (BC); therefore, their influence in the network is not necessarily local but potentially influencing other more distant regions. In conclusion, children with recent onset epilepsy and anxiety demonstrate large scale disruptions in cortical and subcortical brain regions. Network science may not only provide insight into the possible neurobiological correlates of important comorbidities of epilepsy, but also the ways that cortical and subcortical disruption occurs.

Introduction

Anxiety disorders represent a prevalent and problematic interictal psychiatric comorbidity in both adults and children with epilepsy (Beyenburg et al., 2005, Caplan et al., 2005, Ekinci et al., 2009, Goldstein and Harden, 2000, Jones, 2014, Kimiskidis and Valeta, 2012, Kwon and Park, 2014, Reilly et al., 2011, Vazquez and Devinsky, 2003). This relationship has been demonstrated through clinic- (Caplan et al., 2005, Brandt et al., 2010, Ettinger et al., 1998, Johnson et al., 2004, Kwon and Park, 2013), community- (Stefanello et al., 2011) and population-based investigations (Gaitatzis et al., 2004, Kobau et al., 2006, McDermott et al., 1995, Reilly et al., 2015, Tellez-Zenteno et al., 2007). Increased anxiety adversely impacts quality of life (QOL) in epilepsy (Johnson et al., 2004, Choi-Kwon et al., 2003, Kwan et al., 2009), explaining more variance in QOL than traditional clinical seizure variables such as seizure control (Johnson et al., 2004).

The cause(s) underlying elevated symptoms of anxiety and anxiety disorders has been controversial. One view is that anxiety disorders are an expected consequence of the onset, course, and treatment of epilepsy including personal fear of seizures and safety as well as potential unpleasant societal reactions that can result in felt stigma (Asadi-Pooya et al., 2007, de Souza and Salgado, 2006, Gandy et al., 2012, Victoroff et al., 1994). Anxiety symptoms may also be elevated in non-affected persons but who have been exposed to family members with epilepsy such as parents and siblings (Baki et al., 2004, Jones and Reilly, 2016), suggesting that stresses associated with epilepsy may even affect those close to the patient, although the possibility of familial aggregation of anxiety disorders must be considered (c.f., Shimada-Sugimoto et al., 2015 for review).

An alternative view is that anxiety disorders in epilepsy are the consequence of untoward pathophysiological effects of epilepsy on the neural systems that mediate emotional function (Yilmazer-Hanke et al., 2015), a classic example being the relationship between ictally elicited fear states in epilepsy and the amygdala (Cendes et al., 1994, Guimond et al., 2008). Furthermore, both atrophy and hypertrophy of the amygdalae have been reported in adults with epilepsy (Cendes et al., 1994, Minami et al., 2015) and linked to interictal psychopathology (Tebartz van Elst et al., 1999). Increased amygdala volume has also been reported in children with chronic focal (Daley et al., 2008) and generalized epilepsy (Schreibman Cohen et al., 2009) and associated with anxiety disorders and depression.

Also supporting a neurobiological view are findings that behavioral problems in general and anxiety disorders in particular can be detected not only very early in the course of uncomplicated pediatric epilepsies, but may be present antecedent to the first recognized seizure and epilepsy diagnosis (Austin et al., 2001, Jones et al., 2007). Furthermore, recent findings indicate that children with new-onset epilepsy and anxiety disorders exhibit enlarged amygdala (left) and reduced thickness of the left medial orbitofrontal, right frontal pole, and right lateral orbitofrontal regions (Jones et al., 2015), suggesting disrupted circuitry in networks known to mediate anxiety.

While provocative, to date only a very limited, albeit hypothesis-driven set of neural regions have been examined in regard to their association with anxiety disorders. A larger issue is the degree to which diffuse cortical and subcortical networks may be differentially affected in children with epilepsy with anxiety disorders compared to non-anxious children with epilepsy and HC, and here network science may be especially informative in characterizing the complex systems involved (Mears and Pollard, 2016). GT is a mathematical approach to understanding systems as networks. Given that cortical and subcortical structures develop in an organized fashion, covariance analyses between different brain regions can be investigated as a network and analyzed using GT methods, a methodology increasingly used in functional connectivity analyses, white matter connectivity, and covariance analyses of cortical and subcortical volumes in both healthy and medical illness populations (Balardin et al., 2015, Ma et al., 2016, Yeo et al., 2016) including epilepsy (Bernhardt et al., 2016, Gleichgerrcht et al., 2015, Song et al., 2015).

It has been increasingly recognized that brain maturation in childhood involves an organizational process that optimizes network efficiency (Bullmore and Sporns, 2012). The strengthening of network coherence with enhanced cortical thickness covariance in childhood is particularly evident in association cortices compared to primary cortices (e.g. motor, sensorimotor and visual areas), suggesting the orchestration of a more efficient inter-cortical information transfer (Lerch et al., 2006). Epilepsy appears to disrupt this large-scale topology. In children with recent-onset idiopathic epilepsies, cortical volumes covariance was altered with higher network segregation and reduced global integration compared with controls, suggesting alteration of large-scale brain networks. Further, this configuration was more vulnerable to simulated network targeted attacks, implying that this altered network might have fewer parallel or alternative pathways to maintain global integrity (Bonilha et al., 2014). Adults with temporal lobe epilepsy (TLE) also exhibit less network efficiency compared to controls (increased path length and clustering, altered distribution of network hubs) (Bernhardt et al., 2011). It is assumed that anxiety disorder will further compromise brain network configurations in epilepsy. Evidence for this network alteration has only been indirectly inferred from the psychiatric literature. However, there is no direct evidence that anxiety disorder impairs network function in childhood onset epilepsy beyond alterations expected in epilepsy without anxiety.

In this investigation we examined typically developing children as well as youth with recent idiopathic epilepsies with and without anxiety disorders. Network analysis using GT investigated the covariance of diverse cortical and subcortical regions to determine whether and in what way children with epilepsy and anxiety disorders differed compared to non-anxious children with epilepsy and HC. Specifically, we examined the global properties of their networks, the most central/important regions in the configuration of the networks, and the hubs or those areas facilitating communication between pairs of regions. Furthermore, we investigated the network of regions demonstrating high synchronization to amygdala enlargement in the group of children with epilepsy and anxiety in order to determine the influence of those amygdala abnormalities on other neural systems.

Materials and methods

Participants

Study participants included 86 children with recent-onset idiopathic epilepsies and 48 HC aged 8–18 years (Table 1). Inclusion criteria were a diagnosis of epilepsy within the past 12 months, no other developmental disabilities or neurological disorders, normal neurologic examinations, and normal clinical imaging results. A board certified pediatric neurologist (blinded to interview data) confirmed that all participants had focal (idiopathic localization-related epilepsy–ILRE) or generalized (idiopathic generalized epilepsy–IGE) seizures and provided independent confirmation of specific epilepsy syndromes. Participants with focal epilepsies were comprised of children with rolandic epilepsy (22.1%), temporal lobe epilepsy (TLE) (8.1%), frontal epilepsy (9.3%), childhood occipital epilepsy (COE) (1%), and focal epilepsy NOS (10.5%). Participants with generalized epilepsies include juvenile myoclonic epilepsy (JME) (30.2%), absence epilepsy (14.0%), and generalized epilepsy NOS (4.7%).

Table 1

Demographic and clinical characteristics of participants.

	Healthy controls (n = 48)	Epilepsy non-anxious (n = 62)	Epilepsy and anxiety (n = 24)
Age (mean ± SD)	13.30 ± 3.17	12.85 ± 3.50	12.05 ± 2.93
Sex (male/female)	21/27	33/29	9/15
Education	7.02 ± 3.03	6.83 ± 3.65	5.87 ± 2.85
IQ⁎ (mean ± SD)	108.8 ± 10.6	102.1 ± 15.2	101.7 ± 15.6
Syndrome onset (months: mean ± SD)	–	6.63 ± 3.65	7.26 ± 3.89
AED (yes/no)	–	54/8	17/7
Epilepsy syndrome (ILRE/IGE)	–	28/34	16/8

Significant between groups at a p < 0.05. Education represents years of school attendance.

HC were age-matched first-degree cousins of epilepsy participants who presented no history of seizures, current anxiety disorder, any initial precipitating injuries (e.g., febrile convulsions), other developmental or neurologic disease, or loss of consciousness > 5 min. First-degree cousins were used as controls rather than siblings for the following reasons: (i) first-degree cousins are more genetically distant from the participants with epilepsy and thus less predisposed than siblings to shared genetic factors that may contribute to anomalies in brain structure and cognition; (ii) a greater number of first-degree cousins are available than siblings in the target age range; and (iii) the family link was anticipated to facilitate participant recruitment and retention over time. All children and parents participated in a psychiatric diagnostic interview and underwent magnetic resonance imaging (MRI). Further information about participants and the inclusion and exclusion criteria can be found in previous publications (Hermann and Jones, 2006).

Psychiatric diagnostic interview

Every participant and their parents separately participated in a semi-structured interview using the Kiddie–Schedule for Affective Disorders and Schizophrenia–Present and Lifetime Version (K-SADS) (Drevets, 2007). Interviews were recorded and randomly selected for review to ensure inter-rater reliability and a consensus diagnosis was reached when parent and child interviews were conducted by different interviewers. Disorders such as posttraumatic stress disorder (PTSD) and obsessive compulsive disorder (OCD) were not included as anxiety disorders, since recent classification systems (DSM-5 and International Classification of Diseases [ICD-10]) currently place both of these in separate categories. Because of the highly comorbid nature of anxiety disorders, children with other comorbid diagnoses (i.e., attention-deficit/hyperactivity disorder [ADHD] and depression) were not excluded from this study (Kessler et al., 2012). Twenty-four children with epilepsy met the Diagnostic and Statistical Manual of Mental Disorders, fourth edition (DSM-IV) criteria for a current anxiety disorder while the remaining 62 children did not meet such criteria. There were no statistically significant differences in ADHD between groups with epilepsy with and without anxiety as a comorbidity. The anxiety disorders represented included specific phobias, separation anxiety, social phobia, generalized anxiety disorder and anxiety disorders not otherwise specified (NOS). Five of the 48 controls met criteria for anxiety disorder, but given their low number they were not further considered.

Research approval was obtained from the University of Wisconsin School of Medicine and Public Health Sciences Institutional Review Board (IRB). Written informed consent and assent was obtained from parents or legal guardians and participants, respectively, on the day of the study.

Image acquisition and analysis

Images were obtained on a 1.5 T GE Signa MRI scanner (GE Healthcare, Waukesha, WI, U.S.A.). T1-weighted 3D spoiled gradient recall (SPGR) were acquired for each participant (echo time (TE) = 5 msec, repetition time (TR) = 24 msec, flip angle = 40°, slice thickness = 1.5 mm, slices = 124, plane = coronal, field of view (FOV) = 200 mm, matrix = 256 × 256). All MR images are inspected prior to image processing. Image quality (motion, artifacts, signal dropout) is rated on a 5-point scale. We require a grade of 3 or better for a scan to be included in this analysis. Processed data are then visually inspected for proper segmentation boundaries. Editing of the pial and white matter surface is completed when defects are evident.

Images were processed with the FreeSurfer image analysis suite (http://freesurfer.net), a set of software tools for the study of cortical and subcortical anatomy. Volumes of cortical regions were based on the Desikan-Killiany atlas.

Graph theory analysis

Nodes definition and network analysis

Eighty-five nodes were included in the analysis comprising frontal, parietal, temporal, occipital, and subcortical regions. The nodes used and their abbreviations can be found in Supplemental file 1.

Partial correlations between node volumes were calculated controlling for intracranial volume (ICV), which rendered a weighted symmetric 85 × 85 correlation matrix for each group based on the covariance between such volumes. Global metrics were calculated over a range of graph connectivity densities from 5 to 40%, and local measures were calculated at a density level of 26% given that this was the level of full graph connectedness for all groups. The Force Atlas algorithm of the open source software Gephi (https://gephi.org) was used for the 2D visualization of community structure on each group (attraction strength = 10, repulsion strength = 2000, gravity = 30).

Graph theory measures and statistical analyses

GT measures were calculated using the Matlab-based Brain Connectivity Toolbox (BCT, http://www.brain-connectivity-toolbox.net/). In order to investigate statistically significant group differences, each group matrix was resampled by replacement (i.e. bootstrapped) a total of 1000 times. Given that GT measures could be influenced by the nature of the graph (e.g. number of nodes and links), each graph measure was calculated on 1000 random matrices with the same number of nodes and degree distribution as the pertinent graphs. In this way, the null hypothesis could be tested with p-values (corrected for multiple comparisons using Bonferroni correction) for each of the global and local measures. Group differences were investigated by performing ANOVA analyses, and two-sample student's t-test for the global and regional measures, respectively (both correcting for multiple comparisons using Bonferroni correction). No statistical comparisons were undertaken in the visualization of adjacency matrices and community structure since those concerned the actual matrices of the groups.

The characteristic path length, and the clustering coefficient (CC) are well-known GT measures. The former reflects the separation between nodes in a network, and therefore, provides information about network integration. An integrated graph, G is one that allows for information to be efficiently transferred between nodes. However, this measure diverges when nodes are disconnected (have no neighbors), which usually happens at low graph densities. Given that in this work, global GT measures would be acquired over a range of graph densities, the harmonic mean was used instead. The harmonic mean indicates the level of integration of the graph without the problem of divergence; the higher its value the lower the network integration. It is defined as the inverse of the global efficiency, E, which is based on the shortest paths between nodes, N as , where d is the distance from node i to node j (Wang et al., 2010). The CC, which reflects the local segregation in a network, is defined as the ratio of the number of links between each node's neighbors to the total number of links that would exist between them (Watts and Strogatz, 1998). In order to investigate the global segregation, the average CC can be calculated, which is defined as the sum of the local CC of each node divided by the total number of nodes in the network. A network with high level of segregation is one in which intercommunication between nodes' neighbors exist, therefore local transfer of information can be possible. Given that nodes with low degree (low number of connections with other nodes) would have a higher probability for its neighbors to be connected to each other, the average CC would be driven heavily by contributions of lower degree nodes. To avoid this, transitivity, T was used instead, which is the ratio of triangles to triplets in the network: (Newman, 2003, Humphries and Gurney, 2008). The modularity index (MI), informs about the quality of graph subdivision into modules or communities that contribute to the same processes (Blondel et al., 2008, Boccaletti et al., 2006a). In a highly modular network, nodes within the same modules are said to be working toward the same process. In this study, MI was estimated using the modularity Louvain algorithm. These global metrics were calculated over a range of topological thresholds in order to ensure that results are not driven by graph density.

The Z-score with respect to controls was calculated for both epilepsy groups for centrality measures such as strength, eigenvector centrality (EC), and subgraph centrality (SC) at a threshold of 26%. The strength of a node reflects the sum of the weights of the number of connections the node exhibits. A node with high strength is one that has strong correlation inside the network; therefore, having high local influences. EC is a measure indicating the quality of connections in a node. A node with high EC is one that connects mainly to nodes with high degree or high EC themselves (Estrada and Rodriguez-Velazquez, 2005, Spizzirri, 2011); and it is defined as , where a is the adjacency matrix, and λ is a constant. SC evaluates the number of closed walks (loops) in which a node participates while assigning higher weighting to subgraphs formed by fewer nodes (Hermann and Jones, 2006). Nodes with high SC tend to give and receive information in an efficient manner; and it is defined as , where μ is the spectral moment of i. BC is a measure representing the relevance of a node for the communication between other nodes in the network (Boccaletti et al., 2006b). Nodes with high BC facilitate global integrative processes given that they serve as “highways” to facilitate “traffic” flow in the network (Sporns et al., 2007). Therefore, this measure is important when investigating nodes that integrate pairs of nodes that might not be efficiently associated between themselves in other way. It is defined as , where n is the number of shortest paths connecting nodes j and k, and n(i) is the number of shortest paths connecting nodes j and k and passing through i. In this investigation, the network hubs were investigated using BC.

Results

Table 1 provides details of the demographic and clinical characteristics of HC and anxious and non-anxious epilepsy participants. The groups did not differ in age (F(2,131) = 1.164, p = 0.315), sex (χ2 = 2.038, p = 0.361), or grade level (F(2,131) = 0.994, p = 0.373). Both groups of participants with epilepsy presented with full scale IQ that fell in the average range, but significantly lower compared to controls (F(2,131) = 3.759, p = 0.026). Children with versus without anxiety did not differ in age of onset of epilepsy (p = 0.467), number of AEDs (χ2 = 3.178, p = 0.075), proportion with focal versus generalized epilepsy syndromes (χ2 = 3.303, p = 0.074), or rates of comorbid ADHD (χ2 = 0.23, p = 0.43).

Given that total gray matter volume was not significantly different between groups, graphs were calculated while correcting only for ICV.

Matrix visualization

Adjacency matrices showed notable differences in the group with epilepsy and anxiety compared to both controls and non-anxious epilepsy in which subcortical regions were mainly integrated among themselves in the anxious children with epilepsy (Fig. 1). In order to be able to see group differences clearly, we calculated contrast matrices between each group (Fig. 1S in Supplemental document 2), with the least organized differences shown between controls and non-anxious epilepsy.

Fig. 1

Adjacency matrices. Unthresholded matrices for (A) controls, (B) epilepsy without anxiety, and (C) epilepsy with anxiety. The order of nodes is the same as in the Supplemental document 1.

Adjacency matrices. Unthresholded matrices for (A) controls, (B) epilepsy without anxiety, and (C) epilepsy with anxiety. The order of nodes is the same as in the Supplemental document 1.

In order to visualize networks, matrices were thresholded at a density level of 26%, which is the level at which groups showed full graph connectedness. Such 2D matrix visualization can be appreciated in Fig. 2, in which controls show a sparse network when compared to both groups with epilepsy who present highly modular connections. However, the group of epilepsy participants with anxiety show the strongest modular connectivity which is depicted in the separation of same-color nodes. The integration of subcortical structures among themselves in the group with epilepsy and anxiety can also be appreciated in Fig. 2 in which most of the regions that are separated from the rest of the network belong to subcortical areas.

Fig. 2

Community structure. Community structure in controls (top), children with epilepsy without anxiety (middle), and children with epilepsy with anxiety (bottom). The spatial distribution of nodes was calculated using the force-atlas graph algorithm, where nodes that demonstrated stronger connections are located closer in space, while nodes with fewer connections tend to be farther in space. Nodes with a similar color belong to the same module. Calculated at a density of 26%. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)

Global measures

Global measures were calculated from the bootstrapped matrices in each group in order to be able to discern statistical group differences (by performing ANOVA analyses correcting for multiple comparisons). Children with epilepsy and anxiety presented significantly higher transitivity and lower harmonic mean compared to both non-anxious epilepsy and controls (Fig. 3). However, in terms of MI, both groups of children with epilepsy showed similar values that were greater than control participants at each density level. This confirms higher modularity in the children with epilepsy compared to controls while also indicating higher levels of segregation and integration in the children with epilepsy and anxiety compared to the other groups. When testing the null hypothesis on the global measures, each group was statistically significant at each density level (corrected for multiple comparisons).

Fig. 3

Global measures. Transitivity (left), harmonic mean (middle), and modularity index (right). Participants with epilepsy and anxiety presented significantly higher global segregation and integration compared to both healthy controls and participants with epilepsy without anxiety. *Statistically significant between groups after ANOVA testing (corrected for multiple comparisons). Each group result and at each density level is significant (corrected for multiple comparisons) against zero (null hypothesis).

Centrality measures

Measures of centrality were also calculated at a threshold of 26% from the bootstrapped matrices of each group. This allowed calculations of statistical group differences between both epilepsy groups and controls (by performing two-sample student's t-test, correcting for multiple comparisons). These measures represent the most important nodes for the configuration of the network. Fig. 4 shows the ranked Z-score of the node strength in children with epilepsy without (green) and with anxiety (red). As can be seen, non-anxious epilepsy presented higher strength (greater than two standard deviations; bars with different pattern in Fig. 4) than controls in the left medial orbitofrontal gyrus, right fusiform, and left putamen, while the group with epilepsy and anxiety presented high strength in frontal, parietal, cingulate, temporal, and occipital regions compared to controls. Some of the regions presenting lower Z-score than controls were the right pars opercularis and left pars triangularis, bilateral lateral occipital, and left amygdala in non-anxious epilepsy while in anxious epilepsy were the left cuneus, right paracentral, and left lateral occipital gyrus.

Fig. 4

Centrality measures. Ranked Z-score distribution for node strength (top), eigenvector centrality (middle), and subgraph centrality (bottom) on children with epilepsy without anxiety disorder (green), and children with epilepsy with anxiety disorder (red) with respect to controls. Bars with different pattern represent nodes with a Z-score value > 2 or lower than − 2. Calculated at a threshold of 26%. *Statistically significant (corrected for multiple comparisons) against zero (null hypothesis); †Statistically significant between the given group and controls (corrected for multiple comparisons). (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)

Regarding the Z-score distribution of EC, again the left medial orbitofrontal and right fusiform were the nodes showing higher EC in non-anxious epilepsy while again regions in frontal, parietal, cingulate, temporal, and occipital lobes showed high Z-score in the group with epilepsy and anxiety compared to controls (Fig. 4). Again, bilateral lateral occipital, right pars opercularis and left pars triangularis, with the addition of the right paracentral gyrus were the nodes showing lower Z-score for EC in non-anxious epilepsy while for the group with epilepsy and anxiety, again the left cuneus, right paracentral, and left lateral occipital gyrus with the addition of bilateral pars triangularis, and right precuneus presented lower Z-scores for EC compared to controls.

For SC, no nodes showed higher or lower Z-scores compared to controls in non-anxious epilepsy. However, the group with epilepsy and anxiety presented a great number of nodes with high Z-score compared to controls with no nodes showing lower Z-score than controls (Fig. 4).

In summary, the regions presenting the highest Z-score consistently across centrality measures were the left medial orbitofrontal gyrus, and right fusiform for non-anxious epilepsy, and the left rostral middle frontal, left postcentral, left rostral anterior cingulate, left superior temporal, left parahippocampal gyrus, right lingual, and right fusiform for the group with epilepsy and anxiety. Also, those regions consistently presenting lower Z-scores than controls were bilateral lateral occipital, right pars opercularis and left pars triangularis for non-anxious epilepsy while the left cuneus, right paracentral, and left lateral occipital gyrus in epilepsy and anxiety. Since both of these centrality measures convey different information about the node properties and their importance in the network, such consistency proves those nodes to be essential in the configuration of the network at all those different levels of information.

Amygdala connectivity

Given that we previously showed that children with epilepsy and anxiety presented larger amygdala volumes (significantly so for left amygdala) compared to both controls and non-anxious epilepsy (Jones et al., 2015), the network that presented synchronization with bilateral amygdala was investigated in the group of children with epilepsy and anxiety for a threshold of 26% (Fig. 5). Most of the nodes that correlate to amygdala volumes are subcortical structures, with the exception of the left insula and the right frontal pole; both presenting high BC. Given that the left insula and the right frontal pole are part of the nodes with high BC on the network in the group of children with epilepsy and anxiety, their influence in the network is not necessarily local. The other nodes showing high BC in the group with epilepsy and anxiety and on the other two groups can be found in Table 1S in the Supplemental file 2.

Fig. 5

Amygdala network in anxious epilepsy. Amygdala network in children with epilepsy and anxiety. Red = nodes correlated to bilateral amygdala, blue = nodes correlated to the left amygdala, green = nodes correlated to the right amygdala, gray = nodes correlated to amygdala volumes through intermediate nodes. Bigger circles represent nodes with high BC. Calculated at a threshold of 26%. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)

Supplemental analyses

Given that the number of participants in each group differed and there were slight differences in gender, AED, and epilepsy syndrome, we undertook additional analyses to determine if the reported effects were impacted by these factors. In order to do this, the same analyses were repeated using gender, AED, and epilepsy syndrome as nuisance covariates separately. Supplemental analyses (Supplemental file 3) revealed no impact of epilepsy syndrome on the adjacency matrices and the global and local measures, with only minor impact on local measures associated with controlling for AED and gender. Additionally, given that the number of participants in the group with epilepsy and anxiety was lower than both other groups, we also performed the same analyses with only 24 participants in each group (the number of participants in the smallest group—epilepsy with anxiety). Reducing the number of participants per group lowered the power of the analysis and increased the variance, but the same trends prevailed for the adjacency matrices, global measures, and regional measures; the latter presenting minor changes in the nature and ranking of some nodes.

Furthermore, we performed additional analyses based on an atlas with smaller regions (Destrieux atlas from FreeSurfer), in order to confirm that results were not due to the atlas selection. The same trends were found for adjacency matrices, community structure, and global measures. Results can be found in the Supplemental file 4.

Discussion

This study demonstrates the presence of large scale disruptions in cortical and subcortical brain regions associated with the presence of anxiety disorder in children with epilepsy compared to non-anxious children with epilepsy and controls. Furthermore, these disruptions are not limited to anxiety comorbidity since children with epilepsy alone showed significant differences when compared to controls. Since these effects are present near the time of epilepsy onset, the results infer an antecedent disruption of normal neurodevelopmental brain processes in a manner beyond the impact of epilepsy itself as represented by children with epilepsy without anxiety.

The groups of children with epilepsy who were discordant for the presence of anxiety were comparable in intelligence, age of onset of recurrent seizures, number of medications, epilepsy syndrome, as well as the proportion of children with other pertinent comorbidities (e.g., ADHD). Children with epilepsy and anxiety had a current DSM-IV anxiety disorder, the diagnoses of which were determined independently of any knowledge of MRI findings and their analysis. Classification of epilepsy syndrome and determination of other clinical seizure characteristics was made independent of knowledge of the child's psychiatric diagnosis or cognitive assessment.

While there has been controversy regarding the etiology of anxiety disorders in epilepsy, the findings obtained here clearly suggest a powerful effect of disrupted regional brain synchronization linked to the presence of a current anxiety disorder. For the first time, network analysis using GT investigated the covariance of diverse cortical and subcortical regions to determine whether and in what way children with epilepsy and anxiety disorders differed compared to non-anxious epilepsy and HC. The salient findings follow below.

Global network analysis

The network of both groups of children with epilepsy presented higher modularity than controls as shown when the community structure was investigated and when MI was calculated. However, the group with epilepsy and anxiety presented the highest modular organization suggesting non-homogeneous development in epilepsy that is even less homogeneous when anxiety disorder is present. Moreover, modules seem more integrated to each other suggesting high regional interdependence between the different communities. The low modular arrangement observed in controls indicates the opposite. Successful brain development has been linked to widespread cortical thinning (with its associated volume reduction) (Shaw et al., 2006, Shaw et al., 2008), while children with epilepsy have shown a less homogeneous cortical development; therefore, lower modularity in controls compared to patients could be expected. However, similar studies regarding prospective evaluations should be undertaken in order to elucidate the actual developmental trends.

The group of children with epilepsy and anxiety presented both higher segregation and integration compared to both controls and non-anxious epilepsy. This was investigated through examination of harmonic mean and transitivity and could also be inferred via inspection of Fig. 2. Furthermore, adjacency matrices presented similar patterns of connectivity between controls and non-anxious epilepsy compared to anxious epilepsy, the latter showing low synchronization between cortical and subcortical structures like the thalamus, amygdala, pallidum, and caudate. This separation between cortical and subcortical regions can also be observed in the thresholded 2D visualization in which subcortical areas (blue nodes) form a different module.

The right fusiform gyrus consistently exhibited greater Z-score with respect to controls in both groups with epilepsy regarding each of the centrality measures investigated in this study. This was the node presenting high local influences (as indicated by high strength and SC) that at the same time had neighbors with high influences (EC). The fusiform gyrus has been associated with cognitive load in a verbal working memory task in HC (Vogan et al., 2016). However, the fusiform gyrus has also been associated with abnormal resting-state functional connectivity (RS-FC) among anxious subjects compared to HC (Cui et al., 2016). Given that the right fusiform was found to be a highly central region in both groups with epilepsy, its development might be influential over other cortical/subcortical regions, however, it might be differentially affecting each group with epilepsy.

Non-anxious epilepsy also presented the left medial orbitofrontal – a region known to be part of the default mode network (DMN) – consistently exhibiting greater Z-score with respect to controls in strength and eigenvector centrality, while participants with epilepsy showed the same trend for strength and subgraph centrality. The DMN is a network of brain baseline activity that has been proved to be important in self-referential processes, internal mentation, and focused attention (Raichle, 2015); and its disruption has been linked to a diverse number of neurological disorders (Buckner et al., 2008, Greicius et al., 2004, Lustig et al., 2003). Functional and neuroanatomical studies have shown the involvement of DMN-related regions as central measures in healthy brains (Schaefer et al., 2014, van Oort et al., 2014). Given that anxious children with epilepsy did not show high eigenvector centrality of the medial orbitofrontal gyrus compared to controls might demonstrate a potential functional disruption of the DMN in epilepsy when anxiety is a psychiatric comorbidity – a node with high eigenvector centrality is one engaging with other nodes with high quality of connections or potential hubs of the network. Stretton et al. (2015), compared temporal lobe epilepsy (TLE) patients with and without affective disorders, in order to characterize potential differential relationships between DMN activity (task-based functional magnetic resonance imaging) and presence/absence of anxiety. They found a direct relationship between DMN-related regions deactivation and affective disorders, although only 5 out of 17 patients had anxiety as target disorder. The involvement of DMN-related structures in anxiety on healthy individuals has also been suggested (Simpson et al., 2001, Tao et al., 2015), but such relationships have rarely been investigated in adults with epilepsy and never in children with epilepsy.

Contrary to non-anxious epilepsy, the group with anxiety presented a higher number of central regions in frontal, parietal, temporal, and anterior cingulate areas that consistently repeated among centrality measures. The frontal and temporal lobes help mediate language, memory, and executive functions which are often adversely affected in the presence of anxiety disorder (Fuster, 2002). Specifically, working memory (Hopfinger et al., 2000) has been impaired in subjects suffering from anxiety along with sensory-perceptual processes (Cornwell et al., 2007), and processing efficiency (Eysenck and Calvo, 1992). Additionally, the anterior cingulate has been positively associated with cognitive load in HC (Cui et al., 2016) while the primary somatosensory cortex (located in the postcentral gyrus) is important for interoceptive awareness (Khalsa et al., 2009), which is commonly affected in patients suffering from panic disorder (a type of anxiety disorder). Given that these regions were found to be highly central areas in the group of patients with epilepsy and anxiety they might be contributing/spreading neurobiological disturbances given their strategic locations regarding cognition and behavior. The group with epilepsy and anxiety is clearly presenting their most influential nodes throughout the brain while in non-anxious epilepsy are located in both frontal and posterior areas.

Both groups of patients also presented lower Z-scores with respect to controls in lateral occipital regions, however, the group with epilepsy and anxiety also presented the left cuneus as a region with low Z-score compared to controls. This region was found to present lower RS-FC in subjects suffering from anxiety compared to controls (Cui et al., 2016). Therefore, we could speculate that the relatively lower involvement of lateral occipital regions regarding brain volumes synchronization is related to epilepsy, however, such lower involvement of the cuneus might be specific to anxiety disorders. Non-anxious epilepsy also presented lower Z-scores in the right pars opercularis and left pars triangularis while the anxious group also presented consistency across measures in the right paracentral gyrus. Both the pars opercularis and triangularis are well-known regions for language processing in heathy individuals, with their language involvement being found lower in children with epilepsy than in HC (Sepeta et al., 2015). Even though covariance analyses of cortical/subcortical volumes can only inform about indirect influences in regional development, functional associations cannot be discarded since they have been found to be related to brain morphometry (Brandt et al., 2015, Harms et al., 2013).

Amygdala network

When investigating the amygdala network, the results showed that even though a few regions present amygdala synchronization above the selected threshold (Fig. 4), 70% of those regions were hubs in the network (as investigated using BC; see Table 1S of the Supplemental file 2). Given that a hub in the network — as investigated using BC — allows efficient communication between different pairs of nodes, amygdala enlargement is indirectly influencing many brain regions and making its network impact not only within subcortical structures but influencing diverse cortical regions as well.

In this investigation, the right putamen was the only hub common to all three groups which could be the region regulating subcortical structures. Both groups of children with epilepsy shared the left putamen as well as hub, therefore the volume of subcortical areas might be more strictly synchronized in epilepsy. The group with epilepsy and anxiety presented the right frontal pole, the left insula, left rostral middle frontal gyrus, right lingual gyrus, right inferior temporal gyrus, and right pars opercularis as hubs that were specific to that group. The only region where children with epilepsy and anxiety did not present any hubs was in the parietal lobe, contrary to findings in the other groups. Again, the lack of hubs in one of the most important regions of the DMN in children with epilepsy and anxiety might be showing neurobiological disruptions beyond the epilepsy syndrome.

Previous functional neuroimaging studies have demonstrated that patients with social anxiety disorder (SAD) have aberrant functional connectivity, but to our knowledge very few studies have applied GT techniques to this issue in the general population of persons with SAD. Liu et al. (2015) recently examined 20 SAD patients using resting-state fMRI. Whole-brain voxel-wise functional networks were constructed by measuring the temporal correlations of each pair of brain voxels and then hubs were identified by using GT methods. Compared with HC, the SAD patients showed significantly decreased functional connectivity strength (FCS) in bilateral precuneus and significantly increased FCS in the right fusiform gyrus.

In summary, both groups with children with epilepsy presented with abnormal covariance of cortical and subcortical volumes, indicative of antecedent differences in brain development in both groups of patients with new/recent onset epilepsy compared to controls. However, these abnormalities were clearly more exacerbated in a number of ways in the context of anxiety as a psychiatric comorbidity of childhood epilepsy.

Study limitations

It would have been helpful to have a group of normally developing children with anxiety to determine whether the patterns observed in the anxious and non-anxious epilepsy participants were the same or different compared to appropriate controls. While interesting, the primary emphasis here was on the children with epilepsy with versus without anxiety given the existing competing hypotheses regarding the etiology of anxiety in children with epilepsy (social versus biological).

Regarding the methodology used, cortical covariance analyses assume that brain regions that covary together (in thickness, volume, etc.) are developing together. Furthermore, it can be assumed that they could be functionally linked, therefore contributing to similar functional processes. However, this approach can only convey part of the story and does not inform about the direction in which regions positively correlate (e.g. does not tell if regions are increasing or decreasing in volume). Therefore, this can only serve to inform about the positive associations between areas.

Conclusions

Children with recent onset epilepsies demonstrate large scale disruptions in cortical and subcortical brain regions associated with anxiety that go beyond that associated with the presence of epilepsy without anxiety. Network science may serve not only to provide insight into the possible neurobiological correlates of important comorbidities of epilepsy, but also provide insight into the ways that cortical and subcortical disruption occurs. Important are how these networks may further evolve over time, how comorbidity-specific network abnormalities differ, and their eventual resolution following epilepsy remission.

The following are the supplementary data related to this article.

Supplemental file 1:

Nodes

Supplemental file 2:

Additional results

Supplemental file 3:

Supplemental analyses

Supplemental file 4:

Supplemental analyses with a different parcellation atlas

Financial disclosures

The authors report no conflicts of interest.