Global Effects of Focal Brain Tumors on Functional Complexity and Network Robustness: A Prospective Cohort Study

Abstract BACKGROUND Neurosurgical management of brain tumors has entered a paradigm of supramarginal resections that demands thorough understanding of peritumoral functional effects. Historically, the effects of tumors have been believed to be local, and long-range effects have not been considered. OBJECTIVE To test the hypothesis that tumors affect the brain globally, producing long-range gradients in cortical function. METHODS Resting-state functional magnetic resonance imaging (fMRI) data were acquired from 11 participants with glioblastoma and split into discovery and validation datasets in a single-center prospective cohort study. Fractal complexity was computed with a wavelet-based estimator of the Hurst exponent. Distance-related effects of the tumors were tested with a tumor mask-dilation technique and parcellation of the underlying Hurst maps. RESULTS Fractal complexity demonstrates a penumbra of suppression in the peritumoral region. At a global level, as distance from the tumor increases, this initial suppression is balanced by a subsequent overactivity before finally normalizing. These effects were best fit by a quadratic model and were consistent across different network construction pipelines. The Hurst exponent was correlated with graph theory measures of centrality including network robustness, but graph theory measures did not demonstrate distance-dependent effects. CONCLUSION This work provides evidence supporting the theory that focal brain tumors produce long-range gradients in function. Consequently, the effects of focal lesions need to be interpreted in terms of the global changes on functional complexity and network architecture rather than purely in terms of functional localization. Determining whether peritumoral changes represent potential plasticity may facilitate extended resection of tumors without functional cost.

A new philosophy of supra-marginal resection has entered the vocabulary of neuro-oncology whereby surgical excision is extended beyond the confines of a lesion as it is visible on standard structural imaging. This paradigm is based on the ABBREVIATIONS: ANCOVA, analysis of covariance; ANOVA, analysis of variance; ANTs, advanced normalization tools; BOLD, blood oxygenation leveldependent; fMRI, functional magnetic resonance imaging; ME-ICA, multiecho independent component analysis identification of tumor infiltration beyond the contrast-enhancing margin, [1][2][3][4] fluorescenceguided imaging that can identify potentially infiltrated brain at surgery, 5,6 clinical studies correlating extended resections with improved survival, 7-10 and potentially improved effectiveness of adjuvant therapies. 11 However, extending surgery outside the lesion margins comes with the risk of intruding on eloquent brain and consequently the potential of neurological impairment. Therefore, understanding the function of the peritumoral brain is of paramount importance. This requirement to accurately map brain function led to many notable contributions from neurosurgery in understanding functional neuroanatomy. 12,13 However, the traditional approach to understanding the effects of focal tumors per se on brain function is somewhat limited. Firstly, the effects of brain lesions on the surrounding tissue are believed to be only local (ie, through structural disruption), and potential long-range effects are overlooked. Furthermore, function itself is viewed as being uniform and consequently represented as a binary outcome (ie, brain either is or is not eloquent); considering function as a continuous variable has not been investigated. Understanding putative long-range effects of brain tumors and nonbinary gradients in function may allow more accurate prediction of the effects of tumor removal and insight into mechanisms of plasticity-related recovery.
Fractal analysis of resting-state functional magnetic resonance imaging (fMRI) data offers the potential to understanding putative long-range gradients in brain function. Fractals are signals that display scale invariance; that is, they have similar features regardless of the scale at which they are viewed. Fractal signals are pervasive in nature (examples include coastlines, snowflakes, and ocean waves) and in the brain (neuronal membrane potentials, neural field potentials as well as electroencephalography, magnetoencephalography, and blood oxygenation level-dependent [BOLD] contrast endogenous signals). 14 From a fractal perspective, healthy homeostasis emphasizes the complexity of the underlying biological processes rather than regular steady-state behavior. 15,16 In the brain, this complexity is believed to enable adaptation to stimuli or other challenges, while simpler and less complex dynamics are indicative of less advantageous function or disease. [17][18][19][20] Connectomics is a novel multidisciplinary paradigm that naturally views the brain as a complex network of individual components interacting through continuous communication and offers a window into the effects of focal tumors at the global level. 21 Furthermore, graph theory analysis of brain networks unlocks a new terminology with which to describe brain functional neuroanatomy and for predictive modeling of, for example, the effects of tumor removal and how this is related to network robustness. 22 Combining fractal and network analysis methods provides a complimentary approach to identifying putative wide-ranging effects of focal brain tumors on signal complexity and network robustness.
In this work, the concept that tumors only lead to local effects on brain function is challenged, and a new theory that considers the global effects of a tumor is proposed that can augment our current brain mapping techniques. BOLD sensitive, resting-state fMRI was used to map brain function with data-driven analysis of fractal complexity and functional connectivity at the whole brain level. It was hypothesized that tumors would produce long-range effects outside of their margin as depicted on structural magnetic resonance imaging (MRI) and that there would be gradients of BOLD signal and functional connectivity extending away from the tumor. Evidence is presented to support the theory that focal brain tumors produce long-range and nonlinear gradients in function.

Design
The study was approved by the Local Regional Ethics Committee (protocol number10/H0308/23) and was a single-center prospective cohort design. All participants provided written informed consent. Inclusion criteria were the MRI appearance of a tumor consistent with a glioblastoma. All participants had a confirmed glioblastoma at local histological review according to World Health Organization criteria. 23 Resective surgery was performed in all cases under general anesthesia using an operating microscope (OPMI Pentero 900, Carl Zeiss Meditec AG, Oberkochen, Germany) with fluorescence guidance (Gliolan R , medac Pharma, Stirling, United Kingdom) and neuronavigation (StealthStation R S7 R System, Medtronic Inc, Dublin, Ireland).
A cohort of 5 participants was selected to form a discovery dataset for hypothesis generation. Subsequently, a further 6 participants were included, forming the validation dataset. The discovery dataset was designed to be homogeneous and included moderately sized tumors based on the right parietal lobe, whereas the validation dataset was designed to be more heterogeneous with a variety of tumor sizes and locations. Splitting the group into separate discovery and validation cohorts was performed to test if findings were replicable across a variety of tumor presentations. Demographic information is summarized in Table 1.

Imaging
MRI data were acquired using a Siemens Trio 3T scanner and 16channel receive-only head coil (Siemens AG, Erlangen, Germany). A multiecho echo-planar imaging sequence was performed for 10 min and 51 s at a repetition time, TR = 2.42 s resulting in 269, 3-dimensional volumes covering the cerebral cortices and cerebellum. Anatomical images were acquired using a T1-weighted magnetization prepared rapid gradient echo sequence.
A summary of the analysis pipeline is presented in Figure 1. Advanced Normalization Tools (ANTs) cortical thickness pipelines were used for structural analysis (http://stnava.github.io). [25][26][27] A hand drawn tumor mask to encompass the contrast-enhancing margin was created by the neurosurgeon and used as an additional prior for segmentation. Finally, images were nonlinearly mapped to the space of the Montreal Neurological Institute using symmetrical diffeomorphic registration with masking of the tumor from cost-function calculation.

Fractal Analysis
Fractal analysis was performed with a wavelet-based estimator of the Hurst exponent ( Figure 2) applied to processed resting state fMRI data registered to the corresponding structural image. 30 In brief, time series at each voxel underwent wavelet decomposition. A linear model was regressed on to the logarithm of the variance of wavelet coefficients as a function of scale. The estimated slope of this model is proportional to the Hurst exponent, H, and in turn the fractal dimension. Hurst exponents   were obtained for every intracerebral voxel, and then thresholded at H <0.5 (corresponding to white noise) and values were extracted for each parcel of the parcellation template (see below).

Peritumoral Analysis
The hand drawn binary tumor masks were dilated to 30 mm beyond their original margins in 2-mm increments, followed by the subtraction of the preceding masks. This created a series of 2-mm annuli expanding from the tumor margin. These annuli were then masked, only being included in further analyses where they overlapped the cortical grey matter estimated from the structural MRI ( Figure 3). For each mask, for each participant, average H estimated from fMRI was extracted for each annular mask. Masks were subsequently reflected across the interhemispheric fissure to compute the corresponding values from the contralateral, healthy hemisphere. Values were then expressed as the ratio of H in homologous regions.

Connectomics
A summary of the connectome analysis pipeline is presented in Figure 1. Parcellation was initially performed with an anatomical template of 251 equal sized parcels 31 and then subsequently all analyses were repeated on a randomized template of 256 equal sized parcels 32 to test for parcellation independence of network features. Parcels overlapping the tumor mask were removed from subsequent processing. Statistical dependencies between mean processed time series of each parcel were computed using wavelet, Pearson, and partial correlations to test for independence of network features from the underlying method of connectome constructions. Negative correlations were excluded but thresholding was not applied (ie, matrices were 'fully connected').
A graph is composed of nodes (brain regions) and links, which are associated with weights that denote the strength of the connection between nodes. Specific graph theoretical measures were specified prehoc, chosen to reflect fundamental centrality features of the underlying network, including node strength; global efficiency; betweenness centrality 33 ; within-module degree z-score 34 ; participation coefficient 34 ; eigenvector centrality.
A network measure particularly pertinent to neurosurgery is node robustness. This is tested by synthetic lesioning, whereby a node (and all of its connections) is removed from the network and the percentage change (delta) of a selected global network measure is computed. Typically, this is performed with global efficiency in which case the resultant measure is also known as delta efficiency. Each node's delta efficiency values were expressed as a ratio over the value in the homologous node in the contralateral hemisphere. Euclidean distance was defined from the center of each node to the center of the tumor. Network measures were computed using the brain connectivity toolbox and transformed to Z-scores. 35

Statistics
Between-tissue differences in H were tested with Analysis of Variance (ANOVA), with the threshold for significance, P < .05. Testing relationships with distance for the dilating tumor mask and H parcellation analyses were performed with linear (P (k) = a + b * k ), quadratic (P (k) = a + b * k + c * k 2 ), cubic (P (k) = a + b * k + c * k 2 + d * k 3 ), exponential (Pk = a * ek * b ), and exponential including intersection (P (k) = a + b * e k * c ) models. Goodness of fit was determined with the Akaike Information Criteria. 36 Given the distinct recruitment criteria used for the discovery and validation cohorts, differences in the association between log-transformed H and tumor distance where tested using an Analysis of Covariance (ANCOVA) test. All correlations between network measures underwent Bonferroni correction for multiple comparisons.

Tissue Hurst Values
Whole-brain Hurst exponent maps, and H masked by cortical segmentation estimated with ANTs from the structural MRI of the discovery cohort, are presented in Figure 3. Brain extraction, registration and segmentation with ANTs resulted in high-quality images using default parameters without noticeable artifacts due to the presence of the tumor. Differences in the means of the Hurst exponent between neocortical grey matter and white matter for the discovery, validation, and complete cohorts were not significant (ANOVA, F = 1.28, P = 0.31, Figure 4), although the observed higher values in grey, relative to white matter are consistent with the previous findings. 18

Peritumoral Hurst Effects With Distance
Using the tumor mask sequential dilation technique revealed a quadratic increase in H with distance from the tumor (R 2 = 0.98, P < 10 −4 , AIC = −18.2, Figure 5 left) in the discovery cohort. The quadratic increase in H was replicated in both, the validation (R 2 = 0.65, P < 10 −3 , AIC = −17.2; Figure 5 middle) and complete (R 2 = 0.95, P < 10 −4 , AIC = −18.3; Figure 5 right) cohorts. The ANCOVA analysis revealed an interaction effect between log-transformed H and tumor distance (F = 11.7, P < 0.0001) across cohorts. This association was significantly weaker for the validation cohort than for the discovery and complete cohort (P < 0.05; post-hoc Bonferroni corrected). A penumbra of reduced complexity relative to the contralateral hemisphere extended for 10 to 15 mm from the tumor border, while beyond 15 mm complexity continues to increase beyond the values seen contralaterally, albeit nonmonotonically, suggesting a complex relationship to global changes.

Whole Brain Parcellation of Hurst Versus Distance
To expand on the changes in H related to the distance from the tumor border, the H maps were parcellated using the same templates as for the connectome analysis, with the value of each parcel reflecting the mean H. This allowed greater coverage and a whole brain analysis to complement the limited range of the dilating spherical mask, but with coarser resolution than the dilating tumor mask approach. Results were consistent with those from the tumor mask sequential dilation technique in that closest to the tumor values were reduced compared with the contralateral hemisphere, and then increased with distance to above contralateral levels, before finally reducing back to contralateral values ( Figure 6). This trend was best fit with a quadratic model (R 2 0.07-0.12, AIC −9.77 to −12.2; Table 2) in the discovery cohort. The validation cohort also demonstrated this quadratic variation in H with distance from the tumor, again demonstrating replicability as with the sequential mask dilation technique. Furthermore, results were independent of the parcellation technique, with both anatomical and randomized parcellation templates demonstrating the same quadratic variation in H with distance from the tumor in all three cohorts (discovery, validation, and complete). These apparently modest R values reflect that distance is but one of likely many factors driving variation in H. Additionally, the magnitude of these changes is broadly comparable with percentage changes in task-based functional MRI, where signal changes of 2%−4% for primary cortices and 1%−2% for association cortices are typically quoted.

Whole Brain Network Measures Distance Effects
To further expand upon the distance-related changes in H, we investigated whether functional connectome graph theory measures also demonstrated complimentary distance-related effects, specifically that of network robustness. No significant distance-related trends, linear or quadratic, were identified for any graph theory measure (Figure 7). Again, results were consistent across discovery and validation cohorts, anatomical and randomized parcellation templates, and with links weighted by either Wavelet or Pearson product-moment correlations.

Network Measures Relationship to Hurst Exponent
Finally, given that node H, but not graph theory measures were related to distance from the tumor margin, the relationship between H of resting state BOLD time series and functional connectome graph theory measures was explored. H was significantly correlated (r = 0.33-0.87, P < 10 −6 ) with all selected network measures, although distributions were negatively skewed ( Figure 8).

DISCUSSION
In this paper, combined fractal and connectomic analyses were performed in a cohort of participants with brain tumors to determine whether the functional sequelae of focal lesions are apparent at the global rather than just the local level. Fractal scaling showed long-range gradients with distance from the tumor margin that are best fit by a quadratic model that describes proximal suppression of complexity, followed by a subsequent increase, before final normalization to values similar to those in   the contralateral hemisphere. Furthermore, there was a correspondence between fractal properties of the BOLD signal and network measures, specifically those corresponding to centrality and network robustness. Together these findings give credence to the novel hypothesis that the brain responds to focal tumors by reconfiguring its complexity globally, or at least far beyond the tumor margin. Previously, functional localization has been based on the premise that focal injury is responsible for the resultant phenotype; this work suggests that the functional sequelae of focal brain tumors, and the resultant clinical effects, need to be interpreted in terms of spatially extended changes to whole-brain functional complexity and network architecture.
In the context of BOLD signal resting state functional MRI, a higher value of H represents increased complexity of the time series and better neuronal function (ie, a more complex signal is able to carry more information and hence represents increased processing capacity). This has been corroborated in published studies commenting on the relationship between H and cognitive function, neurodegeneration, and pharmaceutical enhancement of performance. [17][18][19]37 In this particular instance, the reduction in the Hurst exponent in the peritumoral region could represent a penumbra of regional cortical suppression, which may arise from direct tumor invasion or distant electrochemical effects. Correlation between centrality measures and network robustness implies the peritumoral suppression of H relates to brain tissue that has already had role reduction within the network. In terms of extended surgical resections, this suggests that removal of the surrounding peritumoral brain might be tolerated at the functional and network levels, and might not result in significant additional functional impairment. Conversely, if the tumor was mediating an electrochemical suppressant affect, focused removal of only the lesion may 'release' this suppression and facilitate functional recovery. Key questions are regarding the clinical effects of this phenomena are how it relates to function and impairment, and whether it recovers after removal of the contrast enhancing tumor mass. In this analysis, no distance-related effects were found for any graph theory measure, but they were in turn significantly correlated with H. This discrepancy in distance-related effects for H and network measures could of course be related to a true lack of correspondence between them, or because the distance-related effects themselves only partially explain the overall variance in the data. H is a local measure of BOLD signal complexity, and this appears to be altered in a way that is consistent across patients. However, these nonlinear distant-related changes to brain physiology may impact graph theoretical measures of centrality that characterize topological properties over spatially extended brain circuits in a more individual manner, diluting consistent distancerelated effects.
Previous applications of resting state fMRI to neuro-oncology have concentrated on mapping primary cortices and determining correspondence with awake brain stimulation findings.
Here we develop these methods in performing whole brain mapping that can find general rules about the brain's reaction to the presence of focal lesions. Task-based functional MRI has developed an established role in clinical practice yet often notes activations beyond that which one would expect a priori. 10 Our data suggest that tumors have divergent effects on the brain depending on distance, implying that multiple mechanisms may be represented by such unexpected activations and not just plasticity.
Effects of focal cerebral lesions on cognition are determined by a number of factors, 38,39 while detailed cognitive examination often reveals highly variable findings between patients despite broadly similar histological and demographic characteristics. 40 Understanding the peritumoral brain at the individual level will be critical to determining the safety of resections, or indeed if it is possible to safely exceed the margins of a lesion in an effort to increase survival. Fractal analysis could be used to identify regions of suppressed cortical function that could be safely resected, or indeed areas that may improve upon resection of the main tumor volume and should be left untouched.
This study is a proof-of-concept investigation into the effects of focal tumors on the complexity of BOLD signals and the functional connectome at the global level. While the overall number of participants was limited, it was felt the advantages of dividing the cohort into discovery and validation groups, namely increased methodological robustness and the ability to test replicability of findings, outweighed any marginal increase in statistical power from analyzing participants in a single cohort. This resulted in relatively small cohorts, but with the benefit of each participant acting as their own control this was not deemed to preclude analysis, particularly given the goal was to identify robust findings at the global level.
High-grade gliomas have previously been noted to have complex effects on the hemodynamic response function and BOLD contrast. 10 We sought to mitigate any potential artifacts with a robust preprocessing strategy including ICAbased denoising. Furthermore, the robust distance-related effects described here suggest an alternative hypothesis whereby the previously described peritumoral effects on BOLD contrast are at least partially related to local cortical suppression rather than hemodynamic or oxygenation abnormalities. Validation of these findings in a cohort of participants with low grade glioma, where the effects on BOLD contrast are not apparent, would further increase the reliability of these findings and is currently being undertaken in a follow-up study. Multimodal imaging (for example, also including FLAIR, Multiparametric Mapping, and NODDI protocols) will be invaluable to fully defining the individual characteristics of specific tissues.
A rational follow-up to this work would be to include cognitive and longitudinal follow-up data to test the generated hypotheses; namely, that cognitive sequelae are dependent on the distance of resection from the tumor margin. Cognitive and longitudinal outcome data will aid inference on the compensatory or decompensatory nature of long-range functional gradients. Finally, larger heterogeneous cohorts will help determine whether these functional effects are consistent or if there are specific drivers underlying individual variation in gradients of changes.

CONCLUSION
Fractal and connectomic analyses suggest global rather than local effects of focal tumors on brain function. Furthermore, exploration of fractals and complex network measures revealed a strong relationship between signal complexity and graph theory measures, notably those reflecting network robustness. It is therefore conceivable that extent of resection, a key prognostic marker for survival, could be approached as an optimization problem balancing increased tumor removal with fractal scaling and network robustness. Further data could be used to test hypotheses that these distance-related effects reflect brain plasticity and potentially allow extended resection of the tumor without functional cost.

Dr Hart is funded by the Wellcome Trust Neuroscience in Psychiatry
Network with additional support from the National Institute for Health Research Cambridge Biomedical Research Centre. The imaging studies were funded by an NIHR Clinician Scientist Fellowship for Dr Price (NIHR/CS/009/011). The views expressed are those of the authors and not necessarily those of the NHS, the NIHR, or the Department of Health. The authors have no personal, financial, or institutional interest in any of the drugs, materials, or devices described in this article.