Connectivity within regions characterizes epilepsy duration and treatment outcome

Abstract Finding clear connectome biomarkers for temporal lobe epilepsy (TLE) patients, in particular at early disease stages, remains a challenge. Currently, the whole‐brain structural connectomes are analyzed based on coarse parcellations (up to 1,000 nodes). However, such global parcellation‐based connectomes may be unsuitable for detecting more localized changes in patients. Here, we use a high‐resolution network (~50,000‐nodes overall) to identify changes at the local level (within brain regions) and test its relation with duration and surgical outcome. Patients with TLE (n = 33) and age‐, sex‐matched healthy subjects (n = 36) underwent high‐resolution (~50,000 nodes) structural network construction based on deterministic tracking of diffusion tensor imaging. Nodes were allocated to 68 cortical regions according to the Desikan–Killany atlas. The connectivity within regions was then used to predict surgical outcome. MRI processing, network reconstruction, and visualization of network changes were integrated into the NICARA (https://nicara.eu). Lower clustering coefficient and higher edge density were found for local connectivity within regions in patients, but were absent for the global network between regions (68 cortical regions). Local connectivity changes, in terms of the number of changed regions and the magnitude of changes, increased with disease duration. Local connectivity yielded a better surgical outcome prediction (Mean value: 95.39% accuracy, 92.76% sensitivity, and 100% specificity) than global connectivity. Connectivity within regions, compared to structural connectivity between brain regions, can be a more efficient biomarker for epilepsy assessment and surgery outcome prediction of medically intractable TLE.


| INTRODUCTION
Brain network disorders lead to changes in structural and functional connectivity between brain regions (Bonilha et al., 2012). These changes are supposed to result from changes within brain regions such as cell death or changes in synaptic connectivity (Jenner & Olanow, 1998;McGlashan & Hoffman, 2000). While cell death, unless being compensated by excess growth of other cell types, can be observed through changes in cortical thickness (Pereira et al., 2012), changes in local connectivity remain elusive. Using high-resolution structural connectivity, where the cortical surface is parcellated into 50,000 regions of interest (ROIs) of comparable size, we test whether brain network disorders are visible at this local level, observing connectivity within brain regions, by looking at the case of epilepsy.
Temporal lobe epilepsy (TLE) as the most common type of epilepsy is characterized by recurrent seizures and can be controlled using antiepileptic medications or surgery (Wiebe, Blume, Girvin, & Eliasziw, 2001). For surgery, locating, and removing epileptogenic zones (EZ) which are involved in seizure generation is crucial; however, seizures can be caused by a complex network rather than a single region (Bartolomei, Chauvel, & Wendling, 2008). Structural connectivity between regions shows distinct changes for patients leading to changes in network dynamics (Hutchings et al., 2015).
However, despite studies on cortical thickness (Pereira et al., 2012), surface area (SA) , or gray/white matter volume (Beheshti et al., 2018;Bernasconi et al., 2004), there have been no studies investigating network changes within brain regions.
In this study, we constructed high-resolution structural networks with around 50,000 nodes that form cortical regions to observe local connectivity within brain regions. We compared our results with a low-resolution global network where 68 cortical regions formed the network nodes. High-resolution networks with 1,000 or more nodes were proposed in the past (Besson et al., 2017;Irimia & Van Horn, 2016; and, using the parcellation into 50,000 nodes on healthy subjects, a modular organization within brain regions could be established . Based on this approach, we investigated how epileptogenic networks change at the local level, within brain regions. Do these local changes also correlate with epilepsy duration (Besson et al., 2017)?
Which cortical regions are most informative when predicting surgery outcome? And, most importantly, is the predictive power for surgery outcome based on high-resolution local network measures more informative than that based on normal low-resolution global networks? 2 | METHODS
Amygdalohippocampectomy and hippocampalsclerosis were confirmed histologically using standard criteria (Kreilkamp, Weber, Richardson, & Keller, 2017). No patients had undergone any previous neurosurgery before mesial TLE surgery. Informed written consent was signed for all participants. Disease duration was computed from age at first seizure onset to age at scan time. For surgery outcome, measured at least 1 year of postsurgical follow-up, patients were classified using the International League Against Epilepsy (ILAE) surgical outcome scale.
2.2 | Data pre-processing and structural network reconstruction T1 images were pre-processed with FreeSurfer (http://surfer.nmr. mgh.harvard.edu). Briefly, the pre-processing steps included intensity normalization, skull stripping, and brain tissue segmentation. Cortical surfaces, such as the pial surface, were extracted and manually checked, edited so that its boundary did not cross the white matter surface boundary. All defects arising from this procedure were manually corrected. We down-sampled the pial surface files with an output ratio of 0.1 using the Matlab toolbox "IsoMesh" (Fang & Boas, 2009) and constructed structural connectomes . The final surfaces were composed of $50,000 triangles, each representing a node in the network and with an average SA of approximately 3.5 mm 2 .
We used "DSI Studio" (http://dsi-studio.labsolver.org) to get deterministic streamline tractography from eddy current-corrected diffusion tensor images which have been processed in FSL (https://fsl. fmrib.ox.ac.uk/fsl/fslwiki/). Generalized q-sample imaging (GQI) method with diffusion sampling length ratio of 1.25, 8-fold orientation distribution function (ODF) tessellation with five peaks resolved on an ODF was chosen when reconstructing FIB files. Directions at corpus callosum were checked to look well (i.e., clean fiber direction at mid corpus callosum and crossing patterns at lateral corpus callosum).
Whole-brain seeding, 0.6*(Otsu's threshold) quantitative anisotropy level, 60 turning angle threshold, 10 mm < tracts length < 300 mm were set and a total of 10,000,000 streamlines was saved. Surface files, streamlines, and "aparc+aseg.nii" files generated by FreeSurfer were linearly registered into the same space. Streamlines whose endpoints terminating within the gray matter were acquired. The high-resolution connectivity matrix was defined by streamline counts between the centers of the closet triangles (Euclidean distance) and was normalized by triangle area. Nodes were sorted according to the default FreeSurfer parcellation, the Desikan-Killany (DK) atlas. Each cortical DK area included hundreds of nodes and can shape an intra- F I G U R E 1 Constructing global (between-area) and local (within-area) structural networks. The final high-resolution connectivity is composed of around 50,000 nodes. Zooming into the high-resolution networks, each within-area structure, for example within the fusiform gyrus (inset), can be considered as a small local network. Note that the local network contained connections among nodes within regions, but contained no connections between regions. The connectivity and 3D visual graphs, generated by NICARA, of an intra-area network, for example, withinfusiform gyrus network are shown here. The low-resolution global connectivity has 68 nodes, each corresponding to a DK atlas area. The 3D subset of low-resolution network connections is also shown by NICARA (in the right of low-resolution network) area network. We counted inter-area streamlines and normalized with a logistic function to create weighted low-resolution networks (68 cortical areas) (Hutchings et al., 2015). Refer to Figure 1 and Supplementary Text S1 in Supplementary Material for more details. To allow for a combination of left and right TLE and investigate ipsi-and contralateral differences, the right regions were flipped in patients who had a right-side surgery (Mueller et al., 2009;Taylor et al., 2018). The demographics of 14 right TLE patients were shown in Table 1.
The same was done with control images. Supplementary Material Table S1 lists the abbreviations of all cortical DK areas. We use the term "area" or "region" to refer to a FreeSurfer (DK) area in the following.
To handle the massive amount of connectivity data that is generated by our within region connectivity workflow, we used NICARA (https://nicara.eu) (Kopetzky & Butz-Ostendorf, 2018)a novel software solution to manage large connectome data sets. We configured the above-described imaging workflow in NICARA so that all external pre-and processing steps could be executed via the NICARA web portal. The resulting brain connectivity was automatically imported, partially analyzed and visualized by the NICARA.

| Network analysis
TLE is often accompanied by surface-based morphometry changes, such as SA reduction in the whole or parts of the brain . Beyond surface abnormalities, structural topological organizations between regions were significantly different in patients.
To investigate and elucidate epileptic brain pathology within and between cortical areas, nine metrics including SA before down-sampling, average fiber length (FL), connectivity strength (S), and six topological network properties (Kaiser, 2011;Rubinov & Sporns, 2010) were considered. For the low-resolution network of connections between 68 regions, we used nodal connectivity strength (S i ), efficiency (E i ), and clustering coefficient (C i ). Average FL was the sum of all of the streamline trajectory lengths divided by the total number of streamlines. Connectivity strength was the sum of connection weights which were normalized by triangle SA for high-resolution networks (Besson, Lopes et al., 2014) and by a logistic function for lowresolution networks (Hutchings et al., 2015). As the high-resolution network is very sparse with about 50,000 nodes and fewer streamlines between nodes, we used binary (unweighted) connectivity to calculate network properties. On contrary, the denser low-resolution network with 68 nodes was weighted in line with previous studies (Bonilha et al., 2015;Hutchings et al., 2015). our sparse high-resolution network, we therefore normalized by network properties of 1,000 randomly rewired networks with the same degree distribution. We also normalized low-resolution network measures by values from 100 rewired networks with the same degree and strength distribution. Finally, we calculated small-worldness (σ). Brain Connectivity Toolbox was used to compute network properties (Rubinov & Sporns, 2010). All measures were corrected for age and gender using a general linear regression model generated from healthy subjects. The relative definitions of all above metrics were listed in Supplementary Material Text S1.

| Within-area local organization analysis
Zooming into high-resolution networks, each DK cortical area can be regarded as a relatively small subnetwork that has 100-1,500 nodes. Regardless of inter-area interactions, the SA of each DK cortical regions was first examined. As the morphometry variation may lead to changes in connectivity, we then regressed out the SA effect by a general linear regression model to pinpoint the influence of other metrics, especially for topological properties (see Supplementary Material Figure S1). Network size (

| Epilepsy duration and surgical outcome analysis
Brain features were proved to change with the duration of epilepsy in patients. Gray matter concentration (gray matter volume divided by total intracranial volume) had a negative association with the duration of epilepsy within the hippocampus, temporal lobe and several limbic structures (Bonilha et al., 2006). Hippocampal volume reduced more if patients had a long history of epilepsy before surgery (Theodore et al., 1999). The alterations of cortical thickness correlation networks intensify over time (Bernhardt et al., 2011). To study epilepsy duration effects, we grouped epilepsy patients into two categories (fewer and more than 20-year duration). At the threshold of 20-year disease duration, two groups had matched group size, gender, surgical outcome, and surgery side distribution, which was beneficial for the following group comparison. For comparison, controls were also categorized into two groups whose age and gender were matched with patient groups. Two control groups were maximum overlapped so that patients can be compared to the same healthy controls. Shorter epilepsy duration is demonstrated to make good surgery outcomes more likely (Bjellvi, Olsson, Malmgren, & Ramsay, 2019). We therefore also examined network differences between patients with surgical outcome. Patients here were classified into two groups: good postsurgical outcome (class ILAE 1-2, 21 patients) and bad (class ILAE 3-5, 12 patients) postsurgical outcome. Refer to Supplementary Material Figure S4, Table S2 for more details.

| Statistical analysis
Network metrics comparison used a two-sided 5,000 permutation test (Zhang et al., 2011) with a statistical significance set as p < 0.05.

| Changes in whole high-and low-resolution networks
In line with previous studies , the whole brain SA was reduced in patients compared to controls (Supplementary Material Figure S5). The majority of cortical regions, in particular around the ipsi-lateral temporal and frontal lobe, had significantly decreased SA in patients. Unlike the low-resolution 68 region networks, the high-resolution networks with around 50,000 nodes showed distinct higher connectivity strength, edge density, and lower global clustering coefficient. Both low-and high-resolution networks observed a decreased average FL (between-nodes) in patients.

| Within-region local changes identification
The changed metrics were plotted in Figure   About half of the affected regions were located in the temporal and frontal lobe. The increase of average FL within-region was dramatic and widespread. However, the average FL across the whole brain was reduced (Supplementary Material Figure S5) which suggests the loss of long-length fibers between regions in patients.
To show examples for within-region local changes, we select two regions which have the most changes in patients compared with controls (cf. Supplementary Material Table S3): the ipsi-lateral PREC and the contra-lateral SMAR (See Figure 3). Compared to controls with a matched number of nodes, the modular organization was reduced with lower clustering coefficient (C) and local efficiency (E local ) in ipsilateral PREC, specifically given that more local neighbor nodes were disconnected (e.g., the network in the left panel of Figure 3(b) is less locally connected with several separate branches and isolated small cliques compared with that in Figure 3(a)). Looking at the connected networks in the brain (i.e., the right panel of Figure 3

| Prediction of surgical outcome
Can the network changes be biomarkers for predicting surgery outcome? For the within-area network shown in Table 2 How does this performance compare to the standard approach of observing changes in the connectivity between regions, that is, our F I G U R E 3 Visualization of local network changes as exemplified for the ipsi-lateral precentral and contra-lateral supra marginal region. From left to right, two subjects (one control, one patient) are shown for both ipsi-lateral precentral gyrus ((a) and (b)) and contra-lateral supra marginal ((c) and (d) accuracy, 80.67% sensitivity and 82.50% specificity for the lowresolution network. The major difference between high-and low-resolution methods was that the outcome prediction, in particular correctly identifying bad surgery outcomes, was much better for the high-resolution within-region network approach.

| DISCUSSION
By using structural connectivity within and between cortical areas, based on DTI, we found significant alterations for TLE that were

F I G U R E 4
The intra-regional change patterns related to duration. Patients were grouped into two parts: one whose epilepsy duration is smaller than 20 years, the others above. For example, the ipsi-lateral entorhinal cortex (ENT), which participates in seizure generation and propagation, showed localized FL increase and decreased connectivity strength. The insular cortex as on the spreading pathways (Isnard et al., 2000) was observed a dramatically increased connectivity strength and higher edge density, which could be the reason for bad surgery outcomes when the insula was not considered as a surgery target (Harroud, Bouthillier, Weil, & Nguyen, 2012;Isnard et al., 2000) . Besides, its higher connectivity strength was a good outcome predictor ($74.05% accuracy in Table 2). The locally changed subtemporal cortex, such as fusiform (FUS, decreased global efficiency) and parahippocampal gyrus (PARH, reduced SA), is also known to be involved in seizures (Alarcon et al., 1997). Changed regions beyond seizure propagation pathways might be indirectly affected by a region that is involved in seizure spreading.

| Changes in TLE and their relation to function
Local network changes might also relate to functional changes in brain connectivity for epilepsy patients. For example, we see network changes in the primary motor cortex, the PREC gyrus, in line with previously reported ictal hyperperfusion due to seizure discharges (Tae et al., 2005). And pathological structural changes occur in PREC gyrus with smaller clustering coefficient and local efficiency indicating destroyed local interaction and lower stability to resist external The SMAR gyrus, a portion of the parietal lobe of the brain, plays a role in phonological processing (i.e., of spoken and written language) and emotional response. Changes in bilateral SMAR gyrus and ipsilateral POPE may be associated with abnormal language processing or a partial displacement of language processing to other regions due to the neuroplastic response to cortical insult, abnormal electrical circuitry, or underlying anatomical anomalies (Devinsky et al., 2000;Hartwigsen et al., 2010).  Note: (I) and (C) represent ipsi-lateral and contra-lateral areas, respectively. 50 repetitions were performed to calculate the mean value and standard deviation. The metrics marked with #(") were significantly smaller (larger) in poor-surgery patients (p < 0.01) compared with good-surgery patients. Others without marks only show the good classification of z-score but a slight difference between two surgical outcome groups. 5,000 permutation tests were used in predictions. The prediction performance of one metric is significant when p < 0.05 and metrics with good predictions were shown in the table.  Table 2 and SA metrics shown in Supplementary Material Table S5 to predict surgical outcome. Table S4). For example, PORB and POPE which were located in inferior frontal lobe and associated with DAN showed increasing changed intensity with duration. Especially, the characteristic path length of contra-lateral PORB was suggested to be a good surgical outcome predictor. Likewise, the INS cortex as parts of SAN, showed large differences in patients and their abnormal intensities get stronger with disease duration. The connectivity strength in the ipsi-lateral INS network was also stronger for poor surgical outcome patients.

| Considerations for surgery
Our analysis of local within-area networks indicated several regions were related to surgery success. Local changes were found to predict epilepsy surgery outcome with an accuracy of about 95.39%. Indeed, the prediction performance is better than the current standard approach of observing global connectivity between regions which has an accuracy of only 88% and better than earlier surgery outcome predictions with 81-88% accuracy (Bonilha et al., 2015;Coan et al., 2016;Sinha et al., 2017). Compared with low-resolution network studies of between-regions, our high-resolution predictive method showed better capability with higher AUC (mean: 0.97 vs. 0.95), accuracy (mean: 95.39% vs. 91.45%), and specificity (mean: 100% vs. 88.67%). The predictive power of poor surgery outcome expressed by specificity in our work is crucial since alternative resection targets could be suggested instead of the planned approach.
While some of the regions related to bad outcome are also regions relate to a longer disease duration, which in turn makes a bad outcome more likely (such as ipsi-lateral INS, contra-lateral LOF, and PORB), other indicated regions are independent of duration. Notably, not all changed regions are facilitating surgery and would therefore be potential surgery targets. Some changes could be due to a relocation of function away from regions affected by seizures that can no longer perform a function. For example, the developing brain can re-establish language organization in the right hemisphere after left-side injury (Krägeloh-Mann, 2004). Indeed, language-related brain regions have a wider spatial distribution in epilepsy patients (Devinsky et al., 2000).  Table 2 and Supplementary

| Brain parcellation and network resolution effect
In our analysis, the DK atlas was used to assign regions to the highresolution network. One of the motivations to select the DK atlas is that it has relatively fewer parcellations (68 regions Figure S2). For an atlas with $1,000 regions, the minimum local network size could be fewer than 10 nodes, which seems not suitable for within-region network analysis. However, given that the DK atlas is based on the brain's geometric information derived from its cortical model, it is unclear if other types of parcellations, such as biological-informed parcellations, would perform better. To explore the effect of the different atlas, the "HCP-MMP" atlas (Glasser et al., 2016) with 360 cortical regions was additionally used to construct a low-resolution network using the NICARA platform. Note that the "HCP-MMP" atlas was created from clustering multi-modal data from the Human Connectome Project (https://humanconnectome. 84.67±3.09%), which indicated a small effect of brain parcellation on the predictive ability of epilepsy surgery outcome by using lowresolution networks.
While the parcellation scheme does not have a big impact on the results, it remains unclear what the effect of node density on the current results as well. Therefore, two additional node densities were tested to construct the structural connectome and re-assess results: one consisted of $25,000 nodes (25 k, average vertex area: $6.3mm 2 ); the other consisted of $12,500 nodes (12 k, average vertex area: $12.6mm 2 ). The structural connectivity comparison at different resolutions was shown in Supplementary Material Figure S6.
Edge density was smaller in the fine-grained network (50 k, that is, and their spatial distribution ( Figure S8(b), (c)). As expected, the regional abnormal patterns for the 25 k resolution network as shown in Figure S8 Table S6), high-resolution networks were better predictors than lowresolution networks and among the three high-resolution networks, while all showing a strong prediction performance, the 50 k-network was the best.
Overall, both parcellation scheme and node density do introduce some variability in our studies. For the low-resolution networks, compared to the structural DK atlas, the multi-modal atlas affects the predictive capacity of epilepsy surgery outcome slightly. While there seems to be no significant difference between the geometric-or biological-informed parcellation for the low-resolution analysis, the DK atlas with larger intra-region networks could be better for analysis of connectivity within-regions. Besides, varying node density from $50,000 nodes ($3.5 mm 2 ) to $12,500 nodes ($12.6 mm 2 ), highresolution networks indicate stronger performance than lowresolution networks. The density of $50,000 nodes which are constrained by the 68-node ROIs (DK atlas) would be preferred as the 50 k-network showed the best prediction and was composed of bigger regions that can benefit other network analyses, such as modularity, in the future.

| Methodological limitations and future studies
Concerning limitations of this study, we created structural connectivity only based on cortical regions and while sub-cortical structures are crucial for seizure propagation, our approach was sufficient to yield biomarkers of epilepsy and of surgery outcome. Secondly, although we found structural network changes, the identification of primary pathology regions should be further examined by other modalities, including EEG and electrocorticography (ECoG) when planning surgical interventions (Lieb, Dasheiff, & Engel, 1991). Thirdly, some region size effects remain even though we eliminated measures that are affected by network size. Although fixing edge density would have been an alternative approach for this challenge in the previous studies (Van Wijk, Stam, & Daffertshofer, 2010), it seems impossible as our high-resolution network is very sparse and would easily lose local structures if connections were removed to yield a fixed edge density or thresholds across subjects. That is also the reason why we only use the binary high-resolution structural connectivity without thresholding to get global network properties. To avoid overfitting for the surgery outcome prediction, we used the k-fold cross-validation method. However, it also has a big risk of overfitting if the group size is small. More patients should be included in the following analysis to examine results further. Evidence has shown the left TLE patients showed larger differences compared to controls than right TLE (Ahmadi et al., 2009;Besson, Dinkelacker, et al., 2014

| Conclusion
In conclusion, TLE patients show characteristic changes in their structural connectivity within regions based on a high-resolution network.
The ipsilateral PREC and the contralateral SMAR gyrus showed large topological changes and some regions, such as the ipsilateral BSTS, PTRI, and contralateral CAC, CMF, LOF, were good predictors of surgery outcome. Specially, ipsilateral INS and contralateral LOF, PORB showed both duration-related changes and good outcome predictive power. Such areas with local connectivity changes might be candidates for surgery in medically intractable epilepsy. Indeed, within-area connectivity was a better predictor of surgery outcome (mean: 95.39% accuracy) than connectivity changes between regions.