Altered functional connectome hierarchy with gene expression signatures in newly-diagnosed focal epilepsy

Objective: Neuroimaging research is providing insights into epilepsy as a disorder of brain connectivity linked to functional impairments which may have an identifiable genetic component. This case-control study aims to identify imbalances in a functional connectome dimension spanning from unimodal to transmodal networks and explore the potential genetic basis of such alterations in patients with newly diagnosed focal epilepsy (NDfE). Methods: We used gradient-based analysis of resting-sate fMRI data comparing cortical gradient maps in patients with NDfE (n = 27) to age and sex-matched controls (n = 36). Using a brain-wide gene expression dataset, gene combinations associated with altered brain regions were then entered into an enrichment analysis. Results: We found an increased differentiation of connectivity profiles between unimodal and transmodal networks in NDfE, which was particularly pronounced in the patients with persistent seizures at 12-months follow-up (n=10). Differences corresponded to gradient score reductions in a visual network and increases in limbic and default mode systems which subserve higher-level cognition. Cortical difference maps were spatially correlated with regional expression of a weighted gene combination. These genes were enriched for disease and ontology terms and pathways previously associated with epilepsy and seizure susceptibility. Interpretations: Large-scale functional hierarchy may be altered from in focal epilepsy from diagnosis and correlate with response to treatment. Combining functional neuroimaging and transcriptional data analysis may provide a framework for understanding the wide-ranging impairments associated with the disorder and mechanistic insight into how gene processes may drive alterations in brain function mediating the genetic risk of epilepsy.


Introduction
Epilepsy is a debilitating neurological condition characterised by unprovoked recurrent seizures effecting around 1% of the population. 1 Developments in both neuroimaging and genetic research has provided insight into underlying pathological mechanisms. 2,3 Despite significant efforts, however, progress in managing epilepsy has been slow, with a third of patients unresponsive to available medication. 4 There is a need for reliable neural markers to optimise care pathways from the onset of the condition, particularly for focal epilepsies which have been most commonly associated with pharmacoresistance. 4 Examining patients with a newly-diagnosed focal epilepsy (NDfE) with advanced neuroimaging methods is still rare, 5 though key to understanding pathophysiological processes independent of effects related to prolonged exposure to anti-seizure medication (ASM) and disease or seizure load.
There is growing consensus that epilepsy is a network disorder, with alterations observed in how discrete brain regions functionally and structurally related to one another. 6,7 Resting-state functional MRI (rs-fMRI) provides a window into spontaneous co-fluctuations of brain activity between regions. Functional network alterations such as increased hubness 8 and local connectivity strength 9 in networks including the deault mode network (DMN) have been identified in focal epilepsy. 10,11 However, despite some promising whole-brain neural circuitry (the connectome) findings, clear and reproducible results remain elusive and do not fully account for the functional impairments identified in focal epilepsy from diagnosis in domains ranging from low-level sensorimotor performance to higher-order cognitive processes. 12,13 Recent work has examined connectomes derived from rs-fMRI with gradient mapping. 14 Gradients represent smooth transitions of connectivity structure along the cortical surface, revealing patterns of large-scale network hierarchy. The principle functional gradient that explains the most variance of brain-wide connectivity captures a transition from unimodal sensorimotor areas to transmodal networks including the DMN associated with higher-order processing. 15,16 Gradient alterations have been identified in a number of clinical populations including autism 17 and generalised epilepsy. 18 Compared to traditional edge-based connectivity, low-dimensional connectome representations show increased phenotypic prediction accuracy across a range of cognitive domains, 19 suggesting potential clinical utility for diagnosing and monitoring neurological conditions. . CC-BY 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.

(which was not certified by peer review)
The copyright holder for this preprint this version posted  https://doi.org/10.1101/2021.07.18.21259977 doi: medRxiv preprint Taking advantage of recent methodelogical advances that allow imaging data to be examined alongside transcriptional datasets, large-scale functional networks have also correlate with variation in cortical gene expression. 20 A study also showed that features of brain structure in patients with psychosis was associated with regional brain expression of schizophrenia-related genes. 21 Given the growing evidence that genetic factors play a key role in the development of focal epilepsies, 22,23 combining functional connectomic with transcriptional data may provide insights into the pathophysiology of epilepsy and clues that lead to novel treatment targets.
The aim of this case-control study was to examine rs-fMRI gradient mapping in conjunction with transcriptional data to: (1) identify alterations in functional network architecture in patients with NDfE; and (2) examine whether such alterations are associated brain expression of epilepsyrelated genes. We expected to find principle gradient score differences between patients (particularly those with poorer seizure outcomes) and controls and that affected regions would be associated with an over-representation of gene patterns/pathways previously associated with epilepsy and seizure susceptibility.

2
Materials and methods

Participants
We analysed data from 27 patients (12 females) with NDfE (mean age, 33.1±11.3, range,  recruited over a 2-year period from clinics at the Walton Centre NHS Foundation Trust, Liverpool, UK. Focal epilepsy was diagnosed by epileptologists based on the assessment of seizure semiology in accordance with International League Against Epilepsy (ILAE) operational classifications. 24 Patients with known progressive neurological disease or provoked, acute symptomatic, unclassified, or primary generalized seizures were not recruited. To increase sample size, we did not limit recruitment to drug-naïve patients. In-depth clinical information on these patients has been previously published. 11 Scanning took place an average of 3.7 months (±2.9, range, 1-11 months) after diagnosis -a time interval when we did not expect adverse medication-related effects on brain function. All patients were followed up one year after scanning to determine current seizure status and response to anti-epileptic medication. All patients were followed up one year after scanning to determine current seizure status and response to anti-epileptic medication treatment. We also scanned 36 neurotypical controls (22 . CC-BY 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity. (which was not certified by peer review) The copyright holder for this preprint this version posted

Imaging analysis
We preprocessed rs-fMRI data using SPM12 ( implemented to remove variance related to confounding signals including movement in the scanner and blood oxygen level dependent (BOLD) response originating in white matter and ventricles using CONN's default denoising pipeline (CompCor). 26 No participants had to be removed from the analysis due to low-quality structural or functional MRI or excessive head motion in rs-fMRI data (translation > 3 mm, rotation > 1°).
. CC-BY 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity. The mean functional time-series from 400 regions, based on the 400 Schaefer parcellation, 27 were extracted for each individual. A 400 x 400 functional brain connectivity matrix was generated for each participant by computing bivariate correlations between each pair of regions with coefficients transformed into z-scores (Fisher z).

2.4
Connectome gradient analysis Fig. 1 shows an outline of the analysis pipeline implemented in this study. We extracted gradients from the functional connectivity matrix of each participant using the BrainSpace Toolbox. 28 The matrix was thresholded with a sparsity parameter of 90 and similarity of connectivity profile between regions was computed with a cosine similarity kernel, following prior approaches. 14 Ten gradient components with decreasing explained variance were computed using diffusion map embedding, a nonlinear dimension reduction approach. 29 In gradient space, regions with similar patterns of functional connectivity are closer to one another than regions with different functional connectivity profiles. The alpha (α) parameter which controls the influence of sampling density was set to 0.5 to retain global gradient space relations between regions. To allow betweenindividual comparisons, Procrustes rotation was used to align gradient components of each participant to a gradient template, computed from the average connectivity matrix of all participants. 30 We focused on the first/principle gradient with the best defined functional implications, as shown in previous work. 17,31 To compare gradient scores mapped back onto the cortical surface in NDfE and control groups, and in PS and SF patient subgroups, we used general linear models (GLM) that controlled for the effects of age and sex. PALM, a non-parametric permutation-based FSL tool, 32 was used for family-wise error rate (FWER) correction of p-values (10,000 permutations).
We also tested for group differences within seven canonical functional networks with PALM (10,000 permutations), correcting for age and sex.

Gene expression analysis
We ran a transcriptomic analysis of case-control gradient differences using the AHBA, a transcriptomic dataset with gene expression measurements collected from six post mortem adult brains aged between 24 and 57 years (human.brain-map.org). 33 AHBA data were collapsed into . CC-BY 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.

(which was not certified by peer review)
The copyright holder for this preprint this version posted July 22, 2021. ; https://doi.org/10.1101/2021.07.18.21259977 doi: medRxiv preprint Schaefer atlas regions and combined across donors using differential stability with the abagen python toolbox (https://github.com/rmarkello/abagen), resulting in a 400 × 15632 regional transcription matrix. Following a previously outlined procedure, 21  and Gene Ontology (GO) enrichment analysis with sub-ontologies (Biological Process -BP, Molecular Function -MF, and Cellular Component -CC). Signed PLS-1 weights were used to rank the input set of genes with Pperm < 0.05. We used default setting which included a p-value cutoff of 0.05 with adjustment, 37 a q-value cutoff of 0.2, a background gene list made up of all genes included in the database, and a minimum and maximum gene set size of 10-500. For more specific insight into synapse function associated with PLS1-derived gene sets, we also explored synaptic enrichment in using SynGO (www.syngoportal.org) with default settings including a background of brain expressed genes.

Gradient results
The first functional gradient explained 14.3% of connectome variance in rs-fMRI data and consisted of a continuum ranging from low-level sensory systems to transmodal networks, with intermediary networks positioned in between (Fig. 2a). Subgroup analyses confirmed that the . CC-BY 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.

(which was not certified by peer review)
The copyright holder for this preprint this version posted July 22, 2021. ; https://doi.org/10.1101/2021.07.18.21259977 doi: medRxiv preprint Connectome alterations in focal epilepsy 8 variance explained was similar in NDfE and HC groups (t = -0.1, p = 0.9). We compared participant gradient scores using surface-based GLMs controlling for age and sex. Globally, we found increased functional separation between unimodal and transmodal networks in NDfE, which was particularly pronounced in the PS group (Fig 2). Differences corresponded to gradient score reductions compared to controls in the visual network, and increases at the DMN end, with decreased density in the middle of the gradient continuum (Fig. 2e). Brain regions with highest principle gradient values in controls tended to show the greatest increases in patients, and conversely, areas with the lowest negative values in healthy controls had the greatest decrease in patients with NDfE. In other words, differences in regional connectivity profiles were maximal at the unimodal and transmodal anchors of the principal gradient, both in the whole NDfE group (r = 0.73, p < 0.001), and in the PS subgroup (r = 0.67 p < 0.001) compared to controls (Fig. 2c).
Group comparisons using non-parametric permutation tests (PALM) revealed increased gradient scores located in prefrontal, temporal and posterior cingulate regions in patient with NDfE compared to HC with reductions in occipital, occipito-temporal and central sensorimotor regions.
Compared to the SF group, scores were similarly increased in a posterior cingulate region in patients with PS compared to those that became SF with reductions located in a visual region. We also summarised surface-wide findings according to the seven functional system decomposition 38 which is part of the Schaefer 400 parcellation scheme (Fig. 2d & 3e). Permutation-based testing (PALM) showed a significant reduction in the visual (t = -2.07, p = 0.023, pCor = 0.12) and sensorimotor (t = -1.78, p=0.042, pCor = 0.224) networks and increases in the limbic network (t = 2.5, p = 0.004, pCor = 0.042) in patient with NDfE compared to HC, though only the limbic network survived FWER correction for comparisons across the seven networks (see pCor).
Gradient score increases were found in the DMN in PS compared to SF patients (t = 2.98, p = 0.002, pCor = 0.019).
2 h e r e

Gene expression related to gradient difference
We used a PLS regression to examine gene expression patterns correlated with cortical maps of case-control differences for NDfE-HC and PS-SF comparisons. The first component (PLS-1) explained 12% of the variance in the data with associated weights positively correlated with cortical NDfE-HC t-statistics (r = 0.58, p < 0.001) (Fig. 3a). In predicting cortical PS-SF t-. CC-BY 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.

(which was not certified by peer review)
The copyright holder for this preprint this version posted July 22, 2021. ; https://doi.org/10.1101/2021.07.18.21259977 doi: medRxiv preprint statistics, the first component explained 11.5% of variance with weights also positively correlated with the predictor (r = 0.41, p < 0.001). We identified 3,258 genes at the chosen permutationderived p-value threshold (Pperm < 0.05). These positive correlations indicate that positively weighted genes are over-expressed in cortical regions where principle gradient scores increased for patients with NDfE and PS, while genes with negative weighting are over-expressed for these groups in regions with decreased gradient scores (Fig. 3c).
3 h e r e

Enrichment analysis of genes related to gradient difference
We found 943 genes with Pperm < 0.05 for the NDfE-HC comparison (263 negatively weighted) (see Dataset1). The protein-protein interaction network for PLS-derived genes from the NDfE-HC comparison had 208 connected proteins and 509 edges/interactions, significantly more than the 334 edges expected by chance (Pperm < 0.001) (Fig. 3d). We also tested the ranked PLSderived gene sets for significant disease, GO and KEGG pathway enrichment (Dataset2). such as hyperalgesia (see Fig. 4). Genes implicated in epilepsy-related disease networks (Fig 5a) were among those that formed the most strongly interconnected cluster of nodes in the PPI network (e.g., KCNC1, SCN1B and GABRA5) (Fig. 3d). The KEGG pathway analysis found significant enrichment in a number of cell signalling pathways including "neuroactive ligandreceptor interaction" (geneRatio = 26/338; p.adjust = 0.048) and "calcium signaling pathway" (geneRatio = 20/338; p.adjust = 0.048), as well as glutamatergic-related mechanisms such as "glutamatergic synapse" (geneRatio = 12/338; p.adjust = 0.048) and "glutamatergic metabolism" is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.
We found 838 genes (282 negatively weighted) for the PS-SF comparison (Dataset1). Notably, no significant disease enriched terms related to epilepsy or seizure were found in this comparison.

Discussion
In this study, we examined whether alterations in large-scale functional networks and the potential genetic basis of such alterations can be identified in patients with NDfE and whether they are associated with future response to treatment (i.e., seizure status) using gradient mapping of rs-fMRI. There is evidence that hierarchy in functional systems such as sensory, motor, and default mode networks is a key principle of brain organisation that determines how information flows across the cortex. 39,40 Gradient mapping linked to what is known about the distribution of gene expression across the brain can offer a more global perspective on network alterations in . CC-BY 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.

(which was not certified by peer review)
The copyright holder for this preprint this version posted July 22, 2021. ; https://doi.org/10.1101/2021.07.18.21259977 doi: medRxiv preprint focal epilepsy with potential to account for the wide-ranging sensory, motor and higher-order cognitive impairments associated with the condition from diagnosis.
Our gradient mapping results showed that the component explaining the most variance in the connectivity data across participants corresponded to a dimension spanning from low-level sensory to transmodal (default mode) systems with other networks in-between [ Fig. 2a], in line with previous work. 39 Crucially, we found increased functional separation between sub-systems in patients with NDfE compared to controls. In our data, this was characterised by lower gradient scores in unimodal regions (visual and rolandic cortices) and higher scores in transmodal hubs (PCC, lateral and medial PFC), with reduced density in mid-gradient networks, such as attentional systems (Fig. 2). These differences were all more pronounced in patients that continued to experience seizures in the following 12 months (PS) after scanning. The general direction of these findings is in line with a recent study that also found principle gradient expansion in patients with genetic generalised epilepsy, 18 in notable contrast to the global contraction of unimodal and transmodal networks observed in patients with ASD. 17 Together, clinical gradient findings provisionally suggest that global expansion may be a broader phenotype of epilepsy or seizure-related brain architecture distinguishable from other brain disorders.
In specific reference to the transmodal end of the principle gradient, there is evidence from rs-fMRI and EEG-fMRI studies that DMN connectivity is altered both in focal epilepsy. 10,41 It is widely acknowledged that DMN core nodes play a key role in internally-focused cognition as well as aspects of language processing and memory, 42 neuropsychological domains also known to be impaired in patients with focal epilepsy. 13,43 In our findings, gradient score difference specifically in the DMN was driven by the PS group (Fig. 2d), indicating that this may represent a promising early marker of pharmacoresistance and treatment outcome. It also leaves open the possibility that response to medication and illness duration could partly contribute to DMN abnormalities identified in previous studies, 44,45 given that these have tended to examine patients with refractory epilepsy. There is also some evidence of dysfunction in networks involved in motor and sensory processes which have also been associated with interictal perceptual as well as cognitive deficits. [46][47][48] Our results suggest that diverse neural and behavioural phenotypes of patients with NDfE could be explained in a more unified framework by alterations in a specific large-scale mode of connectivity, one that determines how information flows from sensory and motor to DMN systems through salience, limbic, frontoparietal and attention-related networks.
. CC-BY 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.

(which was not certified by peer review)
The copyright holder for this preprint this version posted July 22, 2021. ; https://doi.org/10.1101/2021.07.18.21259977 doi: medRxiv preprint Future studies using neuropsychological tests to directly examine the relationship between cognitive performance and functional gradients are advocated and in progress. 49 Intrinsic connectivity gradients have recently been used to capture task-related functional reorganisation of episodic memory processes, as well as language and working memory systems in temporal lobe epilepsy. 50,51 To link our gradient findings of stronger differentiation between unimodal and transmodal systems in NDfE and PS to the histologically measured distribution of gene expression in the brain, we first used a PLS regression to identify weighted gene combinations with cortical expression maps that best matched maps of principal gradient differences between groups (PLS-1). 21 We then ran a series of enrichment analyses to examine diseases and functional processes that have been associated with identified gene combinations in previous genetics and functional genomics research. We found strong correlations between case-control differences of principle gradient scores and cortical gene expression maps for group comparisons (PLS-1). This tight coupling suggests that there is a genetic basis to functional connectome differences between clinical and non-clinical populations and is supported by prior work which found a correlation between large-scale functional networks and variation in cortical gene expression. 20 Comparing patients with NDfE with controls, we found that nearly half of all significant disease enriched terms (8/17) referred to epilepsy or seizure subtypes, indicating that genes previously associated with these terms were over-represented in the weighted gene combinations related to between-group differences in principle gradient scores. Additional investigation revealed that proteins coded by the PLS-derived genes formed a dense PPI cluster, with genes significantly enriched for a number of KEGG pathways and GO processes associated with abnormal neuronal function in epilepsy. For example, we found high enrichment scores in a number of cell signalling pathways known to enhance susceptibility to epileptic seizures such as neuroactive ligand-receptor interaction, 52 glutamatergic processes implicated in aberrant neuronal excitability, 53 pathways involving ion dynamics and features of synaptic transmission, 54 and calcium signaling pathways liked to neuronal synchronization and hyperexcitability. 55 In contrast, no terms specified forms of epilepsy or seizures in the PS-SF group comparison. This makes sense given that this comparison was between patients with gene expression profiles consistent with a diagnosis of epilepsy. Enriched disease terms instead referred to conditions . CC-BY 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.

(which was not certified by peer review)
The copyright holder for this preprint this version posted July 22, 2021. ; https://doi.org/10.1101/2021.07.18.21259977 doi: medRxiv preprint including melatonin deficiency and tumor variants, that together may indicate a higher predisposition for a neuroendocrine/autonomic imbalance with some previous evidence of indirect links to seizures susceptibility. 56,57 A number of membrane transporter processes that mediate drug uptake and efflux were also identified, 58 as well as pathways that regulate intercellular signaling (e.g., Phospholipase D), particularly in cells that are under stress conditions. 59 Significant terms highlighted here could provide a starting point for a targeted, larger scale exploration into response to anti-seizure medication which examines patterns of inter-related diseases, pathways and processes associated with functional and structural brain connectivity hierarchy.
Although the use of a newly diagnosed cohort is a major strength of this study, the sample size was relatively small, particularly for the patient outcome comparison (PS-SF). There are also inherent limitations associated with relying upon gene expression data derived from a different sample (AHBA), though useful brain-specific insights have been gleaned using this approach, 21 which has the advantage of being more easily applied to clinical settings. The replication of findings in a larger sample with associated genome-wide whole blood transcriptome profiling and neuropsychological data will be useful in determining the diagnostic and prognostic potential of combining gradient mapping with gene expression for focal epilepsy.
In conclusion, the findings of this study suggest that large-scale functional networks in rs-fMRI are altered from early stages of focal epilepsy, and correlate with response to treatment. The approach offers a global perspective on network alterations in focal epilepsy with potential to provide a unifying account for the previously characterized wide-ranging sensory, motor and higher-order cognitive impairments associated with the condition from diagnosis. When combined with transcriptional profiling, such alterations can highlight pathways and biological processes relevant to epilepsy and pharmacoresistance. This multi-method approach can contribute to bridging the gap between macro-and micro-scale mechanisms in neuroscience and provide a framework for understanding the interaction between brain processes with the aim of deep phenotyping that leads to precision medicine in the treatment of epilepsy.
. CC-BY 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.

(which was not certified by peer review)
The copyright holder for this preprint this version posted July 22, 2021. ; https://doi.org/10.1101/2021.07.18.21259977 doi: medRxiv preprint Potential conflicts of interest Nothing to report.

Author contributions
CEdB was involved in the conception and design of the study, analysis of data, and manuscript drafting; SSK, AGM, LC and BCB contributed to the design and conception of study and text drafting; BA contributed to acquisition of data.  comparisons of principle gradient scores spatial case-control t statistics were used to predict (e) regional gene expression (Allen Human Brain Atlas) for each ROI in (f) a PLS regression analysis with permutation testing. Genes with Pperm < 0.05 were then examined in (g) functional gene enrichment analyses.   . CC-BY 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity. . CC-BY 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity. (which was not certified by peer review) The copyright holder for this preprint this version posted July 22, 2021. ; https://doi.org/10.1101/2021.07.18.21259977 doi: medRxiv preprint