Cellular Substrates of Functional Network Integration and Memory in Temporal Lobe Epilepsy.

Temporal lobe epilepsy (TLE) patients are at risk of memory deficits, which have been linked to functional network disturbances, particularly of integration of the default mode network (DMN). However, the cellular substrates of functional network integration are unknown. We leverage a unique cross-scale dataset of drug-resistant TLE patients (n = 31), who underwent pseudo resting-state functional magnetic resonance imaging (fMRI), resting-state magnetoencephalography (MEG) and/or neuropsychological testing before neurosurgery. fMRI and MEG underwent atlas-based connectivity analyses. Functional network centrality of the lateral middle temporal gyrus, part of the DMN, was used as a measure of local network integration. Subsequently, non-pathological cortical tissue from this region was used for single cell morphological and electrophysiological patch-clamp analysis, assessing integration in terms of total dendritic length and action potential rise speed. As could be hypothesized, greater network centrality related to better memory performance. Moreover, greater network centrality correlated with more integrative properties at the cellular level across patients. We conclude that individual differences in cognitively relevant functional network integration of a DMN region are mirrored by differences in cellular integrative properties of this region in TLE patients. These findings connect previously separate scales of investigation, increasing translational insight into focal pathology and large-scale network disturbances in TLE.

Introduction

Temporal lobe epilepsy (TLE) is hallmarked by localized pathology of the temporal lobe. It is often accompanied by cognitive disturbances and particularly memory deficits (Meador 2002), which are poorly understood through local pathological markers alone. Indeed, cognition is increasingly thought to depend on the orchestrated functional dynamics taking place on a large-scale structural network of connected brain regions (Stam 2014; Bassett and Sporns 2017). These functional network dynamics can be studied using functional magnetic resonance imaging (fMRI) and magnetoencephalography (MEG), which assess synchronized brain activity of different regions.

The most important properties of brain networks in relation to cognition are segregation and integration (Deco et al. 2015; Cohen and D’Esposito 2016; Horien et al. 2020). Segregation refers to the extent to which nodes are locally interconnected, and integration signifies the level of integrative connectivity taking place, either globally or at a particular node. At specific brain regions, the extent of segregation and integration typically show opposite patterns (van den Heuvel and Sporns 2011; Bertolero et al. 2017): Certain regions have a more segregative topological role (e.g., regions within the visual system), whereas others are considered mainly integrative (e.g., regions in the frontal lobe), indicating that brain regions tend to have a “typical” role in the brain network that can be quantified with either segragative or integrative network measures. For example, centrality indicates the expected level of network integration occurring at any particular node, with nodes showing high integration most likely having low segregation.

Integrative connectivity of the temporal lobe, also overlapping with the default mode network (DMN; Raichle et al. 2001), may be paramount to explain individual memory differences in TLE (DeSalvo et al. 2014; McCormick et al. 2014; Douw et al. 2015). An exemplar study explored the integrative connectivity of the hippocampal circuit with the posterior DMN (Voets et al. 2014) reporting that greater memory deficits were related to lower network integration. An open question remains whether combining different imaging and neurophysiological modalities may improve explanation of cognitive variance in these patients. Multilayer network theory provides a framework to integrate such modalities into a single network consisting of interconnected layers (De Domenico et al. 2013). Previous work suggests that multilayer functional brain network integrative measures of centrality supersede unilayer properties in explaining cognitive functioning in Alzheimer’s disease patients (Yu et al. 2017), but this approach has not been used in TLE or on combined fMRI and MEG networks.

Moreover, the cellular substrates of functional network integration are unknown, limiting our understanding of how focal cellular properties and pathology relate to large-scale network disturbances in TLE and other neurological diseases (Bassett and Bullmore 2009; Stam 2014). “Microstructure-informed connectomics” (Larivière et al. 2019) may connect these scales of measurement (Sporns 2016; van den Heuvel and Yeo 2017). In animals, more integrative structural network regions are comprised of bigger neurons with more axons (Scholtens et al. 2014). Moreover, cross-scale relations between structural brain properties covary with disease characteristics in postmortem studies of multiple sclerosis (Kiljan et al. 2019) and Alzheimer’s disease (Jonkman et al. 2020). However, the cellular substrates of functional network integration as an important correlate of cognitive impairment have been impossible to investigate in vivo.

We leverage a unique cohort of TLE patients undergoing functional neuroimaging and clinical neurophysiological recording as well as tissue extraction of the lateral middle temporal gyrus, a DMN region, through resective neurosurgery (Goriounova et al. 2018). We expected greater network integration to associate with more integrative cellular characteristics in terms of morphology and action potential (AP) kinetics (Poirazi et al. 2003; Eyal et al. 2014; Testa-Silva et al. 2014; Goriounova et al. 2018).

Materials and Methods

Patients

All patients undergoing resective neurosurgery for drug-resistant epilepsy localized in the medial temporal lobe between 2009 and 2016 at Amsterdam University Medical Centers (location VUmc, Amsterdam, the Netherlands) were eligible for participation. All patients underwent temporal lobectomy, which included the lateral middle temporal gyrus. Table 1 gives an overview of clinical information that was extracted from medical chart review. A schematic of the cross-scale analysis pipeline is provided in Figure 1.

Table 1

Patient characteristics

ID	fMRI	MEG	Ephys	Morph	Sex	Age	Lat	Dom	Onset	Dur	Freq	Etiol	Type	AED	IR	DR	WMS
1	A	NA	A	A	F	21	R	Non	11	10	32	LGG	SP	LEV | VPA	49	10	85
2	A	NA	A	A	F	31	L	Non	16	15	13	LGG	CP	CBZ | LEV	55	6	100
3	A	A	A	A	F	29	L	Dom	13	17	32	HS	CP | GTC	LTG | TPM	48	5	NA
4	NA	A	A	A	M	20	R	Non	16	5	2	DNET	CP	CBZ | LEV	52	NA	NA
5	A	NA	A	A	F	27	R	Non	14	13	120	DNET	SP	CBZ | LTG	28	7	90
6	NA	A	A	NA	M	18	R	Non	1	17	150	HS	CP	OXC	NA	NA	NA
7	A	A	NA	A	F	23	L	NA	15	8	9	HS	SP | GTC	LEV | OXC	35	1	82
8	A	A	NA	A	M	53	L	Non	6	47	10	LGG	CP	CBZ | CLB	25	9	65
9	A	NA	NA	A	M	19	R	Non	17	2	450	GLI	CP	CLB | OXC	NA	NA	75
10	NA	A	NA	A	F	35	L	NA	1	34	1	NA	SP | CP | GTC	CZP | LCS | LTG | LEV	36	9	78
11	NA	A	NA	A	M	25	L	NA	9	16	2	HS	SP | CP | GTC	CLB | LCS | LEV	52	12	NA
12	NA	A	A	A	F	31	R	Dom	21	10	4	HS	CP	CBZ | CLB	41	1	80
13	NA	A	A	A	M	49	L	Dom	8	41	2	HS	CP	CBZ | CLB | LEV	45	9	NA
14	A	NA	A	A	F	45	R	Dom	23	22	3	HS	CP	CBZ | CLB | LTG	51	11	101
15	A	NA	A	NA	M	38	R	Non	28	10	6	HS	CP	CBZ	NA	NA	110
16	A	NA	A	NA	M	44	L	NA	35	9	4	HS	SP	LCS | VPA	NA	NA	102
17	NA	A	A	NA	M	53	L	Dom	18	35	4	NA	CP	CBZ	NA	NA	NA
18	NA	A	NA	A	F	30	L	NA	2	28	9	HS	CP | GTC	CLB | OXC	54	1	70
19	NA	A	A	A	M	44	L	NA	4	40	8	HS	CP | GTC	LTG | LEV	30	10	83
20	NA	A	A	NA	F	32	R	NA	6	26	151	MCD	NA	NA	NA	NA	NA
21	A	A	A	NA	F	41	R	Non	8	33	6	GLI	AB | CP	LEV | CBZ | VPA | CLB	44	NA	86
22	A	NA	A	NA	M	29	R	Non	23	7	5	NA	AB | CP	PHB | CBZ	28	NA	NA
23	A	A	A	NA	F	20	R	Non	4	17	6	MCD	CP	LEV | CBZ	58	NA	NA
24	NA	A	A	NA	M	21	L	NA	8	13	6	GG	AB	CBZ | LEV | LTG	NA	NA	84
25	NA	A	A	NA	F	48	R	NA	14	34	48	HS	AB | GTC	ZNS | CBZ | VPA	48	NA	104
26	A	NA	A	NA	M	44	R	Non	25	19	1	HS	CP	LTG | CBZ	NA	NA	NA
27	NA	A	A	NA	M	43	R	NA	6	37	7	HS	SP | CP	OXC | LEV	NA	NA	NA
28	A	A	A	NA	F	33	R	Non	19	14	150	HS	AB | CP	CBZ | LEV | CLB	51	NA	80
29	A	A	A	NA	M	51	R	Non	4	48	60	NA	CP | GTC	CBZ | PHB	41	NA	NA
30	A	A	A	A	F	18	R	Non	5	13	17	GG	CP	OXC	66	NA	NA
31	A	A	A	A	M	23	R	Non	11	12	8	HS	CP	OXC	56	NA	83

Abbreviations: ID = identification number, Ephys = electrophysiology, Morph = morphology, Lat = lateralization of epileptic focus, Dom = dominant hemisphere of resection, Onset = age in years at seizure onset, Dur = duration of epilepsy in years, Freq = frequency of seizures per month, Etiol = etiology of epilepsy, Type = type of epilepsy, AED = anti-epileptic drugs, IR = immediate recall on the Rey Auditory Verbal Learning Test, DR = delayed recall on the Rey Auditory Verbal Learning Test, WMS = Wechsler Memory Scale Verbal Memory Index, A = available, NA = not available, F = female, M = male, R = right, L = left, Non = non-dominant, Dom = dominant, LGG = low-grade glioma, HS = hippocampal sclerosis, DNET = dysembrionic neuroepithelial tumor, CAV = cavernoma, GLI = gliosis, MCD = malformation of cortical development, GG = ganglioglioma, SP = simpel partial, CP = complex partial, GTC = generalized tonic clonic, AB = absence, LEV = levetiracetam, VPA = valproic acid, CBZ = carbamazepine, CLB = clobazam, OXC = oxcarbazepine, LTG = lamotrigine, LCS = lacosamide, PHB = phenobarbital.

Figure 1

Schematic representation of multi-scale analyses. In (A), cellular tissue collection from the middle temporal gyrus (pink node) is depicted. Morphological analysis and electrophysiological recordings were performed. In (B), the functional network measure of eigenvector centrality is illustrated for the unilayer (fMRI or MEG) and multilayer (combined fMRI and MEG) network analyses.

All procedures were performed with the approval of the Medical Ethical Committee of VUmc, and in accordance with Dutch license procedures and the Declaration of Helsinki. Written informed consent was provided by all subjects for data and tissue use for scientific research.

Data Availability Statement

Patients did not consent to sharing their raw data. However, the full derived variable set used for analysis in the current work is available from GitHub (https://github.com/multinetlab-amsterdam/projects/tree/master/multiscale_integration), where our code to do so can also be found.

Memory Functioning

Patients underwent cognitive assessments during presurgical workup, as previously reported (Goriounova et al. 2018). In order to assess verbal memory functioning, the Wechsler Memory Scale (WMS) and the Dutch version of Rey’s Auditory Verbal Learning Test (RAVLT) were selected for analysis. From the WMS, the composite Verbal Memory Index was used, which has a mean value of 100 with a standard deviation of 15 in a normative and healthy population. From the RAVLT, immediate recall in terms of total correctly encoded words across five trials (15 words each) was quantified, with possible scores ranging between 0 and 75 words in total. In addition, delayed recall of the 15 words was assessed after 15 min, yielding an additional outcome measure of memory retrieval. Higher values indicate better memory performance.

Structural Neuroimaging

MRI was performed on a 1.5 T magnet (Siemens Sonata) and included an anatomical 3D T1-weighted MPRAGE scan (sequence parameters: time repetition [TR] = 2700 ms, time echo [TE] = 5.2 ms, time to inversion [TI] = 950 ms, 1-mm isotropic resolution, 176 slices). Image processing was performed using FSL5. Standard procedures were used to preprocess structural imaging: Non-brain tissue was removed from the 3D T1-weighted images using the Brain Extraction Tool, and gray and white matter segmentation was performed using FAST. Non-brain tissue was removed and tissue segmentation was performed. To construct each individual’s functional brain network, the Automated Anatomical Labeling atlas was used to define 78 cortical regions. This atlas was warped from standard space to native space, and masked with each subject’s native gray matter mask.

Functional Magnetic Resonance Imaging

Patients underwent presurgical fMRI for language localization. Previous studies have shown that the intraindividual effect of any task or state on network topology is small in comparison with interindividual differences in network topology (Gratton et al. 2018; Kraus et al. 2021), indicating that we may use these pseudo resting-state data to investigate individual differences in network integration. Resting-state analysis has previously been used on such task data for network analysis in both healthy controls and patients (Harris et al. 2014; Krienen et al. 2014; Derks et al. 2017).

Scanning was performed using a standard echo-planar imaging sequence (TR = 2850 ms, TE = 60 ms, 144 volumes, 3.3-mm isotropic resolution, 7 min). During the scan, patients performed a language task in which nine volumes of word generation were alternated with nine volumes of rest (imagery of a landscape).

Preprocessing was performed using standard procedures (Beckmann et al. 2005), including discarding the first five volumes, motion correction, spatial smoothing, and high-pass filtering. Six regions with low signal quality (mainly orbitofrontal areas) were excluded, leaving 72 cortical regions for analysis. In addition, ICA-AROMA was applied to minimize the impact of movement (Pruim et al. 2015). Functional images were co-registered to the anatomical scans. Time series were extracted from the centroids of all regions, after which a 72 × 72 connectivity matrix per subject was created using Pearson correlation coefficients. Finally, the absolute values of these correlations were used as weighted connectivity.

Magnetoencephalography

Patients underwent resting-state MEG as part of their presurgical work-up and/or in the setting of scientific research (Nissen et al. 2017, 2018). Interictal eyes-closed recordings were acquired in supine position using a whole-head system (Elekta Neuromag Oy) with 306 channels inside a magnetically shielded room (Vacuumschmelze GmbH). Data were recorded with a sampling frequency of 1250 Hz, filtered online with a 410 Hz anti-aliasing filter and a 0.1-Hz high-pass filter. The head position relative to the sensors was recorded continuously with head-localization coils. A 3D digitizer (Fastrak, Polhemus) digitized the head-localization coil positions and scalp outline (roughly 500 points) to allow surface matching of the scalp surface points with anatomical MRI.

Three eyes-closed resting-state recordings of typically 15 min each were recorded for clinical analysis of interictal epileptiform activity. Only one recording was analyzed in this study and chosen according to the following criteria with descending priority: 1) consisting of at least 5 min of data, 2) displaying the smallest number of artifacts as per visual inspection, and 3) being the earlier dataset of the three recordings.

Further analysis of these data has been extensively described before (Nissen et al. 2017). Offline spatial filtering of the raw data removed artifacts using the temporal extension of Signal Space Separation (tSSS) using MaxFilter software (Elekta Neuromag Oy; version 2.1). The reconstruction of neuronal sources was performed with an atlas-based beamforming approach, after which time series (virtual electrodes) for each centroid of each atlas region were reconstructed (Hillebrand et al. 2016). These time series were then filtered in the theta band (4–8 Hz), because of its proven relevance for cognitive functioning in these patients by others and ourselves (van Dellen et al. 2009; Douw et al. 2010; Jin and Chung 2015), and to limit the number of investigated variables in this limited sample.

As a measure of functional connectivity, the phase lag index (PLI) was used. The PLI assesses the phase-relationship between two regions by quantifying the asymmetry in the distribution of instantaneous phase differences between two time series. It is robust against zero-lag phase synchronization due to volume conduction or field spread. This analysis yielded a 78 × 78 MEG connectivity matrix per patient.

Network Analysis

We then created minimum spanning trees (MSTs) to extract the functional backbone of each fMRI and MEG network and adequately allow for comparison between patients and modalities without having to use a subjective threshold for pairwise connectivities (Stam et al. 2014). The MST is a binary network that connects all nodes in a network without forming cycles while maximizing connectivity strength. Eigenvector centrality was then calculated per brain region. Of note, there is a plethora of network measures that measure integration, which are highly intercorrelated (Oldham et al. 2019). Instead of using multiple measures in this sample with limited statistical power, we chose to focus on eigenvector centrality as the network measure of integration. Eigenvector centrality is a spectral centrality measure, which not only takes into account the number of connections of a node, but also weighs the number of connections of its neighboring nodes (Lohmann et al. 2010).

In addition to modality-specific analysis of centrality, we also explored multimodal centrality, since an open question is whether combining different imaging and neurophysiological modalities may improve explanation of cognitive variance in these patients. Multilayer network theory offers an analytical framework that allows for such synergy between modalities to be captured (Mucha et al. 2010). Each layer in a multilayer network represents a network characterized by one type of connectivity. Interlayer connections link the same region across different layers. Importantly, multilayer network measures supersede summed properties of individual layers when trying to explain the behavior of other types of complex networks (Stegehuis et al. 2016). Multilayer techniques have only recently been applied to neuroscience, but show promising results towards explaining more cognitive variance than unilayer analyses in for instance Alzheimer’s disease patients (Yu et al. 2017). We therefore used the fMRI and MEG MSTs to construct a two-layer network consisting of the 72 identical nodes representing each brain region in the atlas available in both modalities. Each region was connected only to itself across the two layers, forming an interconnected binary multiplex network. We then calculated multilayer eigenvector centrality per region (De Domenico et al. 2015). Ultimately, this network analysis yielded three centrality values per patient (fMRI, MEG, and multilayer).

All analyses were performed using in-house Python scripts (publicly available from our GitHub page (https://github.com/multinetlab-amsterdam/projects/tree/master/multiscale_integration) in combination with the publicly available Brain Connectivity Toolbox (Rubinov and Sporns 2010) implemented in Matlab R2012a (Mathworks).

Single Cell Electrophysiology and Morphology

Tissue exclusively originated from the lateral middle temporal gyrus and was removed in order to gain access to the disease focus in deeper lying structures. The tissue was resected prior to temporal lobectomy, from the middle temporal gyrus at the location 4-cm posterior to the temporal pole as a block of ~ 1–1.5 cm in diameter. The location of resection was consistent across patients independent on the extent of temporal lobectomy, and can be found in Figure 1 of Goriounova et al. (2018). In all patients, the resected neocortical tissue was not part of the epileptic focus or tumor and displayed no structural/functional abnormalities according to presurgical MRI investigation and histological analysis by an experienced pathologist. Data analysis has been extensively described before (Goriounova et al. 2018).

Upon surgical resection, the cortical tissue was immediately transferred to ice-cold artificial cerebral spinal fluid, then transported to the electrophysiology lab within 15 min, where neocortical slices (350-μm thickness) were prepared (Goriounova et al. 2018). Whole-cell patch-clamp recordings were made of layers 2 and 3 pyramidal neurons and APs were elicited by incrementing step current injections (step size 30–50 pA). Waveforms were sorted according to their instantaneous firing frequency (1/time to previous AP) and AP rise speed was defined as the peak of AP derivative (dV/dt) for all APs from all neurons from each subject.

During electrophysiological recordings, cells were loaded with biocytin through the recording pipette. After recording, slices were fixed in 4% paraformaldehyde and cells were revealed with chromogen 3,3-diaminobenzidine tetrahydrochloride using the avidin–biotin–peroxidase method. Neurons were digitally reconstructed using Neurolucida software (Microbrightfield). Only neurons with virtually complete dendritic structures were included.

We then selected three representative properties pertaining to integration at the cellular scale. Larger dendrites may enable neurons to have more synaptic contacts, putatively playing a more important integrative role than neurons with smaller dendrites (Poirazi et al. 2003; Eyal et al. 2014). Larger dendrites also directly influence the speed of AP initiation, possibly providing these cells with better temporal resolution and more efficient information transfer (Eyal et al. 2014; Testa-Silva et al. 2014; Goriounova et al. 2018). Thus, regions that act as integrators for cognitive processes may be characterized by neurons with larger dendrites and faster APs. We therefore extracted total dendritic length (TDL) of all basal and apical dendrites and then averaged these data (1–10 neurons per patient) as our first measure of cellular integration, also because this measure proved cognitive relevant in these patients before (Goriounova et al. 2018). In addition, we selected rise speed of the first AP, and APs fired at frequencies between 20 and 40 Hz for cross-scale analysis, also due to their proven relevance for cognition in this patient cohort (Goriounova et al. 2018).

Statistical Analysis

Statistical analysis was performed in Matlab R2012a (Mathworks).

Pairwise associations between functional network centralities and memory were tested using Spearman’s correlation coefficients, because of the small sample size and likely non-normal distribution of our data. For these associations between memory functioning and network centrality, we corrected for the nine tests performed across different measures.

Cross-scale pairwise associations were also tested using non-parametric Spearman’s correlation coefficients because of the sample characteristics, with bootstrapping (1000 samples, 95% confidence intervals [CI]). When significant, robustness was explored by permuting the micro–macro pairs to create a data-specific correlation distribution (10 000 permutations) and by leave-one-out analysis. Spatial specificity of the associations was explored by correlating the cellular measure with functional network centralities of all other ipsilateral nodes in the network. The threshold for statistical significance was set at two-tailed alpha < 0.05, but we also report significant results after applying Bonferroni correction for multiple comparisons. In the cross-scale analysis, we corrected for the nine associations tested between network centralities and cellular properties. For the spatial specificity analysis, we corrected for the 36 comparisons made between cellular properties and regional network centralities.

Results

In 15 of the 46 TLE patients originally included (Goriounova et al. 2018), MEG and fMRI were not available. Therefore, 31 patients (15 females) with a mean age of 33 years (± 11 years) were included for the current analysis. Table 1 describes age, sex, clinical information on histopathology, disease duration, side of seizure focus, etcetera, per patient.

Functional Network Integration and Memory Functioning

First, we sought to confirm that functional network centrality of the lateral middle temporal gyrus, as part of the DMN, related to memory functioning in these patients. Greater fMRI network centrality was significantly related to greater RAVLT delayed recall (rho = 0.857, CI [0.373 1.000], P = 0.024, n = 7; Fig. 2), although this correlation did not survive correction for the nine tests performed. Unfortunately, very few patients had complete data across these scales, rendering the remaining pairwise testing of associations underpowered.

Figure 2

All pairwise correlations. This figure shows an overview of all associations between cellular properties and functional network centrality (top three rows), as well as between functional network centrality and memory functioning (bottom three rows), using the maximum samples of patients with available data. In addition, the middle three rows (in gray) reflect all pairwise correlations between scale-specific properties when considered in the same subgroup with complete functional network data for the morphological (n = 5) and electrophysiological (n = 7) analyses. Green elements reflect positive correlations (with rho, P and n in text), ranging from small correlations (rho < 0.4, light green) to medium correlations (0.4 < rho < 0.6, medium green) to large correlations (rho > 0.6, dark green). White, bold text indicates statistical significance (P < 0.05) and very dark green elements with white, bold text indicates statistical significance after Bonferroni correction for multiple (nine) comparisons. Blue elements reflect negative correlations in the same way. RAVLT immediate = Rey Auditory Verbal Learning Test immediate recall, RAVLT immediate = Rey Auditory Verbal Learning Test delayed recall.

Cross-Scale Correlations of Integrative Properties

We then asked whether greater functional network centrality of the resected region went hand in hand with longer TDLs and faster APs, hypothetically signs of greater integrative potential at the cellular level (Fig. 2). Longer TDL was significantly related to greater fMRI network centrality (rho = 0.758, CI [0.115 0.975], P = 0.016, n = 10; Fig. 3). This association remained significant when creating a sample-specific distribution of the correlation through permutation analysis (rho cutoff = 0.636, Fig. 3 left panel). We then performed leave-one-out analyses, where we iteratively excluded a single patient from the correlation analysis to see whether this result was driven by individual data points. This yielded significant results in 7 of 10 analyses (Fig. 3 middle panel). Finally, we explored whether TDL of the resected region also correlated with network centrality of other cortical regions. This analysis revealed a significant negative correlation between TDL of the resected region and fMRI network centrality of the ipsilateral paracentral lobule (rho = −0.721, P = 0.024, Fig. 3 right panel), although this result did not survive correction for the 36 correlations tested. Regional TDL therefore seems specifically correlated to fMRI network centrality of the same region.

Figure 3

Significant pairwise associations between cellular properties and functional network centrality. (A) Displays a scatter plot of the significantly positive correlation between TDL (in mm) and fMRI centrality. The dark blue diamond in the plot represents a single patient, who has low centrality and short dendritic length (schematically depicted in the middle panel), whereas the light blue diamond represents a patient with high centrality and long dendritic length (right panel). In (B), the analyses of the robustness of these results are displayed. The left panel shows the distribution of permuted correlations (10 000 permutations), with the actual (pink) association being smaller than the alpha = 0.05 threshold indicated by the dotted line. The middle panel displays the 10 leave-one-out associations in addition to the real correlation, yielding a significant correlation in 7 of 10 analyses as indicated by the asterisks. The right panel displays all correlations between ipsilateral functional network centrality values (n = 36 regions) and TDL of the resected area. The reported positive association (with the resected region) in pink as well as the negative correlation between paracentral centrality and TDL are significant. In (C), a scatter plot of the significantly positive correlation between AP rise speed (first) and fMRI centrality is shown. In (D), the analyses of the robustness of these results are displayed in the same manner as in (B).

Furthermore, faster AP rise speed (first) was significantly related to greater fMRI network centrality (rho = 0.539, CI [0.042 0.868], P = 0.041, n = 15; Fig. 3 and ). This finding remained in permutation testing (rho cutoff = 0.525), but was mostly nonsignificant in leave-one-out analyses. Spatially, AP rise speed (first) of the resected area was also significantly associated with fMRI network centrality of the precentral region (rho = 0.582, P = 0.025). Moreover, there were significant negative relationships with fMRI network centrality of the superior parietal region (rho = −0.657, P = 0.010) and temporal pole (rho = −0.611, P = 0.018), although these correlations did not survive correction for multiple comparisons and may thus indicate regional specificity of this association.

Although the multilayer analyses were considered exploratory due to the small sample size, the only pairwise correlation that survived correction for the nine pairwise tests performed between cellular and functional network characteristics of the resected region was the significant positive association between AP rise speed (first) and multilayer network centrality (rho = 0.964, CI [0.698 1.000], P < 0.001, n = 7; Fig. 4). Moreover, this association was robust in permutation and leave-one-out analyses, and was spatially specific, even when not correcting for multiple comparisons.

Figure 4

Significant pairwise association between cellular properties and functional network centrality after correction for multiple comparisons. (A) Displays a scatter plot of the significantly positive correlation between AP rise speed (first) and multilayer centrality, which is the only pairwise association significant after Bonferroni correction for multiple comparisons. The dark blue diamond in the plot represents a single patient, who has low centrality and slow APs (schematically depicted in the middle panel), whereas the light blue diamond represents a patient with high centrality and fast APs (right panel). In (B), the analyses of the robustness of these results are displayed. The left panel shows the distribution of permuted correlations (10 000 permutations), with the actual (pink) association being smaller than the alpha = 0.05 threshold indicated by the dotted line. The middle panel displays the seven leave-one-out associations in addition to the real correlation, yielding a significant correlation in all analyses as indicated by the asterisks. The right panel displays all correlations between ipsilateral network centrality values (n = 36 regions) and AP rise speed (first) of the resected area. The reported positive association (with the resected region) in pink is the only significant correlation.

Since different numbers of patients were available for each pairwise correlation depending on modalities available, we also investigated whether the small subset of patients with complete data showed the same patterns of correlation as the larger but heterogeneous samples. Indeed, Figure 2 also shows the associations of both unilayer and multilayer network indices with cellular measures in the same samples of five patients (for TDL) and seven patients (for AP rise speeds). As was the case for the mixed samples, only multilayer network centrality and AP rise speed (first) were significantly correlated in this complete dataset, indicating that patient selection differences between the above described subsamples did not confound the reported pairwise correlations.

Post Hoc Analyses

In order to explore whether our main finding of a positive correlation between functional network centrality and cellular integrative properties held up across different network measures, the Brainnetome atlas (BNA; Fan et al. 2016), and MEG frequency bands, post hoc analyses were performed, the details of which can be found in the Supplementary Materials.

In order to interpret this large number of additional analyses, we here focused on those associations with a medium to large effect (rho > 0.4). As can be seen in Supplementary Figure 1, fMRI degree centrality revealed similar positive associations with cellular properties that we found for eigenvector centrality across the two atlases. Significant was the positive correlation between BNA fMRI degree centrality and rise speed of the first AP (rho = 0.642, P = 0.010). MEG theta band degree centrality was also positively associated with cellular properties when using the AAL atlas, significantly so for TDL (rho = 0.677, P = 0.016), but not for the BNA atlas.

The clustering coefficient yielded mostly negative correlations with cellular properties, as could be expected, revealing a significant negative correlation between BNA fMRI clustering coefficient and TDL (rho = −0.673, P = 0.039).

Eigenvector centrality based on the BNA atlas yielded only small (rho < 0.4) associations with cellular properties, for both fMRI and MEG networks, whereas BNA multilayer eigenvector centrality did yield medium positive associations (0.4 < rho <  0.6) with cellular properties. Analysis of eigenvector centrality in the additional MEG frequency bands revealed several medium positive correlations (0.4 < rho <  0.6) with cellular characteristics.

Discussion

We report on a unique cohort of TLE patients with data spanning multiple scales of investigation, in order to explore the cellular substrates of cognitively relevant functional network integration of a DMN region. As expected (McCormick et al. 2014; Voets et al. 2014; Douw et al. 2015), greater functional network integration of the lateral middle temporal gyrus related to better memory functioning in these patients. Moreover, greater network centrality of this brain region correlated with signs of greater cellular integration: Patients with longer dendrites and faster APs also displayed more integrative functional network profiles.

Our findings signify that brain organization in terms of integrative propensity is preserved across scales of measurement in TLE patients. Functional network centrality of the investigated DMN region correlated with neuronal morphology in terms of TDL of its constituent neurons. In previous work, we showed that pyramidal neurons of patients with higher intelligence in the related cohort had larger, more complex dendrites (Goriounova et al. 2018). Larger dendrites would enable neurons to have more synaptic contacts, putatively playing a more important integrative role than neurons with smaller dendrites (Poirazi et al. 2003; Eyal et al. 2014). Greater dendritic length and axonal density of neurons at particular locations have previously been linked to structural network integration as measured with diffusion MRI: more integratively connected network regions tended to have bigger neurons with more axons, whereas locally connected network regions were made up of smaller neurons that were connected with lower axonal density (Scholtens et al. 2014; Kiljan et al. 2019). The current results suggest that functional network integration also mirrors more integrative neuronal morphology, at least in the lateral middle temporal gyrus as part of the DMN.

Functional network centrality also related to AP kinetics: more integrative functional network topology related to faster AP rise speeds. We previously found that patients with higher intelligence had faster APs during high frequency activity, in addition to the abovementioned longer TDLs (Goriounova et al. 2018). Of note, larger dendrites also directly influence the speed of AP initiation, putatively offering more integrative power (Eyal et al. 2014; Testa-Silva et al. 2014; Goriounova et al. 2018). Our current results put these local indicators of integration into a cross-scale network perspective, as we report on associations with functional network centrality. Our findings support a scale-free view on these type of brain properties, such that within-region cellular integrative properties are reflected by between-region functional network integration, which is of particular relevance towards understanding (cognitive deficits in) TLE. Based on our current findings, one may hypothesize network centrality of the DMN is a mediator between cellular properties and memory functioning. Future studies with larger sample sizes may further elucidate whether this is indeed the case.

Of note, the functional network correlates of local AP kinetics became particularly evident when taking both functional modalities (fMRI and MEG) into account through multilayer network analysis. This finding suggests that multimodal functional centrality may capture cellular brain properties better than either of the two modalities alone, which would be in line with other multilayer brain network studies (Battiston et al. 2017; De Domenico 2017). It should, however, be noted that the lack of robust significance for unilayer network centralities as well as the significance of multilayer network centrality may reflect limited statistical power in small subsamples. Although it is reassuring that all centralities (unilayer and multilayer) show the same positive association with TDL and AP rise speed, it should be noted that particularly the multilayer analysis sample was small (n = 7) and results may not generalize robustly.

Two tests of memory functioning were explored: the verbal memory index of the WMS (Wechsler 2009) and the RAVLT (van der Elst et al. 2005). Although the WMS did not show any significant correlations in this study, the RAVLT delayed recall did. There are both differences and similarities between these tests: They both measure verbal memory encoding, although one is focused on simpler list learning (RAVLT) and the other on learning word pairs (WMS). Previous work in healthy subjects suggests that scores on these two types of tests are only moderately intercorrelated (Wechsler 2009; Miller et al. 2012; Thiruselvam et al. 2015), possibly due to the generally free-recall format of the RAVLT, with potentially richer learning occurring in the WMS due to the feedback that is given after each trial. Future work should address the specificity of our results for these different types of verbal memory.

Based on our findings, richer hypotheses may be formulated on the relationships between local cellular organization and functionality on the one hand, and large-scale network topology on the other hand in TLE. Next steps may include investigation of cellular and network properties of pathological brain regions in addition to the non-pathological area of the DMN assessed in this study. This would bridge the gap between focal TLE pathology and its widespread cross-scale effects on the brain, ultimately pertaining to cognitive impairment. Such studies are of particular interest, since our current finding of cross-scale preservation of integrative propensity in a non-pathological brain region may not necessarily be specific to TLE: Previous work in macaques and postmortem donors also report on such correlations (Scholtens et al. 2014; Kiljan et al. 2019; Jonkman et al. 2020), raising the question whether cross-scale integration is a basic organizational principle conserved across species to begin with. By also involving focal pathological cellular properties (and their large-scale network counterparts), disease-specific processes may be disentangled from such fundamental principles, ultimately leading to a more comprehensive framework on the development of cognitive complaints in TLE. Such a framework could contribute to the formulation of new treatment targets for the cognitive impairments that many patients suffer from. For instance, neurostimulation to remaining parts of the DMN may induce cellular and functional network changes that go hand in hand with cognitive improvement (Meisenhelter and Jobst 2018).

Several limitations of this study should be noted. First and foremost, the unique nature of this multi-scale dataset meant that only a small sample was available, also precluding subgroup analysis or exploration of confounders that might have affected the current results, for example, lateralization of the epileptogenic zone, sex, and age. It was also impossible to use a control group for this analysis: Healthy individuals can obviously not be subjected to the cellular measurements we report on. Also, although resting-state data was available for MEG, we used pseudo resting-state fMRI data, which may have induced deviations from previous works. However, interindividual differences in functional network topology, which were our focus, are several orders of magnitude larger than task versus rest differences within individuals (Gratton et al. 2018). Taken together with the fact that we found positive associations between centrality and cellular integrative measures for both fMRI and MEG, it is unlikely that task performance during fMRI impacted our main findings, but future studies may detail the way in which network connectivity during task performance relates to both resting-state connectivity and cellular properties and memory functioning in these patients. Another limitation is the spatial resolution of matching between the cellular and functional network analyses: Although certainty about the location of origin of the resected tissue was in the order of millimeters, the atlas region used to reflect the tissue location spans several centimeters. We chose this atlas as it has been used successfully by our group in previous studies with comparable patient populations and thus has proven cognitive and cellular relevance (van Dellen et al. 2012, 2014; Carbo et al. 2017), but future studies may aim to match tissue locations with increased spatial accuracy. Furthermore, in post hoc analyses, we confirmed the positive association between fMRI and multilayer (degree and/or eigenvector) centrality and cellular properties when using the BNA with more and functionally defined parcels, but could not replicate these same associations for MEG centrality alone. These findings could be due to the intrinsically lower resolution of MEG, as it is estimated based on adaptive parcellation that no more than ~70 cortical parcels capture the MEG connectome (Farahibozorg et al. 2018). As such, the centroid voxel of the BNA region at the tissue location may have been less representative of the more widespread source of the MEG signals in this area. As for MEG frequency bands, we focused on the theta band due to previously reported associations with memory functioning in these patients (van Dellen et al. 2009; Douw et al. 2010; Jin and Chung 2015), although other frequency bands may also relate to memory functioning. In our post hoc analyses, we also found several non-significant but medium positive associations between centrality and cellular properties in the alpha and broadband frequencies. Future studies may further explore how exactly the different frequency bands relate to cellular properties.

In conclusion, we show that individual differences in functional network integration of a DMN region relate to cellular morphology and AP kinetics. These results underline the translational nature of individual differences in brain properties between TLE patients, which has clinical relevance in terms of memory functioning. Ultimately, such “microstructure-informed connectomics” (Larivière et al. 2019) may lead the way towards better understanding and treatment of neurological disease in general, and memory functioning in TLE patients specifically.

Click here for additional data file.