White matter network topology relates to cognitive flexibility and cumulative neurological risk in adult survivors of pediatric brain tumors.

Adult survivors of pediatric brain tumors exhibit deficits in executive functioning. Given that brain tumors and medical treatments for brain tumors result in disruptions to white matter, a network analysis was used to explore the topological properties of white matter networks. This study used diffusion tensor imaging and deterministic tractography in 38 adult survivors of pediatric brain tumors (mean age in years = 23.11 (SD = 4.96), 54% female, mean years post diagnosis = 14.09 (SD = 6.19)) and 38 healthy peers matched by age, gender, handedness, and socioeconomic status. Nodes were defined using the Automated Anatomical Labeling (AAL) parcellation scheme, and edges were defined as the mean fractional anisotropy of streamlines that connected each node pair. Global efficiency and average clustering coefficient were reduced in survivors compared to healthy peers with preferential impact to hub regions. Global efficiency mediated differences in cognitive flexibility between survivors and healthy peers, as well as the relationship between cumulative neurological risk and cognitive flexibility. These results suggest that adult survivors of pediatric brain tumors, on average one and a half decades post brain tumor diagnosis and treatment, exhibit altered white matter topology in the form of suboptimal integration and segregation of large scale networks, and that disrupted topology may underlie executive functioning impairments. Network based studies provided important topographic insights on network organization in long-term survivors of pediatric brain tumor.

Introduction

Cancers of the brain and central nervous system are the second most prevalent type of cancer in children (Ostrom et al., 2017). Medical and technological advances in cancer treatments have resulted in improved survival rates for children with brain tumors, and research has shifted to emphasize quality of survival, identify psychosocial and neurobiological factors that predict poor outcomes, identify protective factors that promote resilience, and develop effective interventions to address the problems that arise as survivors age and reach adulthood (Lassaletta et al., 2015; Mulhern et al., 2004; Murdaugh et al., 2017).

Research on long-term outcomes of adult survivors of pediatric brain tumors has demonstrated significant detrimental cognitive effects especially in the domain of executive functioning (Edelstein et al., 2011; Hocking et al., 2015; McCurdy et al., 2016a; McCurdy et al., 2016b; Spiegler et al., 2004; Wolfe et al., 2012). Studies have also found that deficits in executive functioning underlie deficits in social and adaptive functioning (King et al., 2015a; Wolfe et al., 2013). As a potentially modifiable domain with consequences to adaptive skills and functional outcomes, executive functioning is an important area of inquiry in the study of long term outcomes in adult survivors of pediatric brain tumors.

Neuroimaging studies in survivors of pediatric brain tumor have found reduced white matter integrity in specific tracts (Aukema et al., 2009; King et al., 2015b; Mabbott et al., 2006; Palmer et al., 2012; Riggs et al., 2014; Rueckriegel et al., 2010; Smith et al., 2014) and that the integrity of these tracts is associated with performance on neurocognitive measures of processing speed, motor speed, full scale IQ, fine motor function, working memory, memory, and executive functioning (Aukema et al., 2009; Brinkman et al., 2012; Law et al., 2011; Palmer et al., 2012; Riggs et al., 2014; Rueckriegel et al., 2015). These studies suggest that structural changes to white matter may represent the neurobiological underpinnings of outcomes. These structural neuroimaging studies in survivors so far have used a univariate framework to identify regions of the brain that differ between survivors and healthy peers, based on the assumption that discrete regions of the brain are responsible for specific functions.

The brain, however, functions as interconnected networks of neuronal information from distributed areas. Knowledge on such complex networks cannot be achieved solely by studying individual components of the networks, as was done in prior literature. Frameworks that emphasize connectivity between brain regions can provide complementary information to traditional neuroimaging techniques, especially when examining higher-order behaviors such as executive functions that depend on the integration of information from spatially distributed regions of the brain. Complex network analysis, such as graph theory, emphasizes how each brain region is connected to others as a system and how ensembles of brain regions work together in a unified network (Rubinov and Sporns, 2010; Sporns, 2012). In graph theory, graphs are composed of nodes and edges; when applied to analyze neuroimaging data, each node is defined as a specific region of the brain, while edges are defined as the existence or strength of the connections between each node to every other node. Graph theory is then used to obtain quantitative metrics describing system level properties. These metrics can be divided into measures of segregation (i.e., the extent to which information is processed locally within a small region), measures of integration (i.e., the extent to which information is processed across spatially distributed regions), and measures of nodal importance via centrality (i.e., properties of nodes that describes its importance within the network based on the spatial locations of the nodes) or vulnerability (i.e., properties of nodes that describes its importance within the network based on perturbations).

Research studies of clinical populations of patients with stroke, schizophrenia, Alzheimer's Disease, traumatic brain injury (TBI), epilepsy, and multiple sclerosis have primarily found alterations in measures of segregation, integration, centrality, and disruptions in small-world properties in chronic stages of injury (Aerts et al., 2016; Bullmore and Sporns, 2009; Crofts et al., 2011; Reijmer et al., 2015; Tang et al., 2015). Moreover, clinical presentations and the degree of executive functioning impairment in various disorders are related to the integrity of brain networks (Bullmore and Sporns, 2009; Caeyenberghs et al., 2014; Panigrahy et al., 2015; Yuan et al., 2015). Disruptions to hub regions also are an important feature of many clinical conditions. Hubs, which have exceptionally high connections to other nodes in the network, are a feature of healthy human brains and are essential for integrative processing (Hagmann et al., 2008; Li et al., 2013; van den Heuvel and Sporns, 2011; van Straaten and Stam, 2013). Disruption of hub regions have been demonstrated in network studies of clinical populations and may represent a final common pathway in the disease process of all neurological disorders (van den Heuvel and Sporns, 2013). In chronic stages of injury, these regions often show significant decreases in measures of centrality across a variety of different brain disorders (Crossley et al., 2014). Further, the level of disruption in hub regions are related to behavioral outcomes (Caeyenberghs et al., 2012; Fagerholm et al., 2015; Kim et al., 2014; Yuan et al., 2015). Taken together, these results suggest that graph theory metrics are sensitive to structural changes that occur as a result of neurological insult, and that they demonstrate concurrent validity with behavioral measures.

Graph theory methods have yet to be used to examine white matter network properties of survivors of childhood brain tumors but are well suited for exploring this population, as there are a multitude of disease and treatment factors that result in white matter disruption (de Ruiter et al., 2013; Edelstein et al., 2011; Gragert and Ris, 2011; Ris and Noll, 1994). For instance, radiation treatment results in dose-dependent and progressive white matter damage (Connor et al., 2016). Chemotherapy is associated with reductions in white matter volume (Ren et al., 2017). Hydrocephalus, a common neurological condition in brain tumor survivors, results in increased intraventricular pressure, damage to periventricular white matter and possible axonal degeneration (Krishnamurthy and Li, 2014). Surgical resection of the tumor beyond the tumor margin as often required causes a loss of brain tissue in areas of the brain that are in the same neural pathway but distal to the site of the lesion (Ailion et al., 2016). Tractography studies in adult patients have found displacement and tumor infiltration in white matter tracts as a result of fast-growing tumors (Nilsson et al., 2008; Wei et al., 2014) with persistent effects of thinning, interruptions and reductions in the tract size after surgery (Lazar et al., 2005). As these treatment and brain tumor related factors contribute to white matter disruption, and since these white matter disruptions are hypothesized to underpin cognitive and functional impairments, a network analysis framework is particularly apt when studying adult survivors of childhood brain tumors.

Two aims were proposed to examine white matter network properties of adult survivors of pediatric brain tumors. The first aim was to establish whether white matter topology is altered in adult survivors compared to healthy controls. It was hypothesized that measures of integration and measures of segregation would be lower in survivors when compared to controls. It was also hypothesized that hub regions would exhibit reductions in measures of centrality and would be preferentially impacted. Finally, it was hypothesized that risk factors such as younger age at diagnosis, longer time since diagnosis, and higher levels of neurological/treatment risk would be associated with more changes to measures of integration, segregation and hub centrality. These variables have been identified as risk factors for worse cognitive outcomes in the domains of intelligence, executive functions, attention, and working memory (Briere et al., 2008; Edelstein et al., 2011; Spiegler et al., 2004). Notably, these risk factors are also related to lower overall white matter integrity and lower white matter integrity in specific tracts (King et al., 2015b; Law et al., 2011; Reddick et al., 2014). The second aim was to establish whether white matter network topology related to behavioral performance on executive functioning (specifically cognitive flexibility) measures. Research supports that executive functioning relies on frontal-subcortical systems, rather than any one region (Bonelli and Cummings, 2007; Koziol and Budding, 2009; Koziol and Lutz, 2013; Mega and Cummings, 1994). Its reliance on the integrity of the system makes using graph theory approaches particularly relevant when relating to behavior. It was hypothesized that lower levels of integration, segregation, and higher levels of overall hub disruption would be correlated with worse cognitive flexibility. It was also hypothesized that differences in cognitive flexibility between survivors and controls would be mediated by global network properties. Finally, it was hypothesized that the relationship between cumulative neurological risk factors and cognitive flexibility would be mediated by topological properties of the white matter network.

Methods

Parent study and procedures

Participants for this study were recruited as a part of a parent study investigating long-term functional outcomes in survivors. The parent study was reviewed and approved by the local institutional review board, and all participants provided written informed consent. Participants were recruited through opt-in letters and newsletters and screened over the phone to ensure they were over the age of 17 and at least 4.5 years post diagnosis. Participants were excluded if English was not their first language, if they met criteria for pervasive developmental disorders, if they were diagnosed with neurofibromatosis, or if they had experienced any other significant neurological insult unrelated to the brain tumor. Eighty-eight survivors met initial criteria for the study and were invited to participate.

On the first day, participants were interviewed to obtain medical history. Written permission was obtained to access medical records to corroborate diagnosis and treatment. Participants underwent a comprehensive neuropsychological evaluation and were provided breaks throughout the day to minimize fatigue. Finally, they were screened to determine whether they could safely participate in the neuroimaging part of the study on a different day. Fifty-one individuals participated in the imaging part of the study, while the other 37 survivors either could not participate due to MRI safety exclusions, indicated that they were not interested, or were lost to follow-up. Manual inspection of the raw diffusion data indicated that three survivors exhibited excessive (i.e., visible and prolonged) head motion across imaging slices and ten had imaging artifacts from surgical implants that led to high levels of signal loss and image distortion. The remaining 38 survivors with good quality of diffusion imaging data represented the sample for this study.

Healthy adults were recruited to serve as the comparison group for analyses through the research pool of Department of Psychology at Georgia State University as well as fliers and advertisements in the community. Participants completed screening for MRI safety and were administered the SCID-II (First et al., 1997) to ensure that they did not currently meet criteria for psychological or substance abuse disorders. A total of 58 healthy controls were scanned, of which three were excluded due to severe head motion. Of the remaining 55 control participants, 38 participants were selected to comprise the control sample with age, gender and handedness matched with those in the survivor group. All controls had no history of a neurological illness. Control participants were administered the same comprehensive neuropsychological battery on the first visit and completed the one-hour imaging portion of the study on a different day following the same procedure used in controls. The average interval between the evaluation and the MRI scan for the entire sample was 90 days. The average interval for survivors was 86 days (SD = 195), while the average interval for the controls was 93 days (SD = 179). The two groups did not significantly differ in the average amount of time from the evaluation to the MRI scan, t(74) = 0.16, p = .88.

Survivors were paid $100 for the time and travel associated with partaking in the neuropsychological and imaging part of the study. Community participants were paid the same amount, while participants recruited from the psychology department pool received class credit for the neuropsychological testing part of the study and $50 for the imaging part of the study.

Participants

Demographic and treatment characteristics of the samples are described in Table 1. Mean age, gender, and socioeconomic status were not significantly different between the two groups (p > .05). The control group had higher levels of education, higher IQs, and were more ethnically diverse than the survivor group.

Table 1

Demographic, diagnostic and treatment characteristics.

	Controlsn = 38	Survivorsn = 38	Group Difference Statistic
Demographic information
Females (n, %)	21 (55%)	21 (55%)	χ2(1, N = 76) = 0, p = 1.0
Ethnicity (n, %)			χ2(2, N = 76) = 10.18, p < .01⁎
Caucasian	13 (34%)	29 (76%)
African-American	14 (37%)	4 (11%)
Latino/a	4 (11%)	2 (5%)
Asian	5 (13%)	1 (3%)
Mixed	2 (5%)	2 (5%)
Socioeconomic statusa			χ2 (1, N = 76) = 3.45, p = .06
High	21 (55%)	28 (74%)
Middle/low	17 (45%)	9 (24%)
Mean age at examination (SD)	22.54 (4.83)	23.11 (4.96)	t(74) = 0.50, p = .62
Mean years of education (SD)	14.47 (1.98)	13.39 (2.39)	t(74) = 2.14, p = .04
IQ scaled score (SD)	111 (9)	98 (18)	t(74) = 3.33, p < .01
Diagnostic information
Mean age at diagnosis (SD)		9.03 (5.03)
Mean years since diagnosis (SD)		14.09 (6.19)
Range (years)		4.5 to 30
Tumor type (n, %)
Medulloblastoma		12 (32%)
Low-grade astrocytoma		13 (34%)
High-grade astrocytoma		1 (3%)
Craniopharyngioma		2 (5%)
Ganglioglioma		3 (8%)
Ependymoma		2 (5%)
Otherb		5 (13%)
Tumor location (n, %)
Posterior fossa		26 (68%)
Temporal lobe		4 (11%)
Occipital lobe		1 (3%)
Fronto-parietal lobe		2 (5%)
Temporal-parietal lobe		1 (3%)
Hypothalamus		1 (3%)
Medulla		1 (3%)
Third ventricle/sellar/suprasellar		2 (5%)
Treatment information
Hydrocephalus (n, %)		25 (66%)
Radiation treatment (n, %)		20 (53%)
Chemotherapy (n, %)		15 (40%)
Endocrine disorder (n, %)		20 (53%)
Neurosurgery (n, %)		37 (97%)
Total resection		26 (68%)
Subtotal resection		11 (29%)
Seizure medications		3 (8%)

Note. Intelligence was measured by the first version of the Wechsler Abbreviated Scale of Intelligence (Wechsler, 1999).

Given the small cell sizes in some of the ethnicity categories, the chi square test was carried out using three levels for the Ethnicity variable: Caucasian, African-American and Other (the combined participants in the Latino/a, Asian, and Mixed categories).

SES = Current socioeconomic status was calculated using the Hollingshead Four factor Index of Social Status (Hollingshead, 1975). High SES consisted of scores 1 and 2 on the scale, while Middle/Low SES consisted of scores 3, 4, and 5 on the scale.

1 Oligodendroglioma, 1 choroid plexus papilloma, 2 PNET Not Otherwise Specified, 1 Mixed astrocytoma/ganglioglioma.

Measures

Measure of cognitive flexibility

The Color Word Interference Test from the DKEFS was used as the measure of cognitive flexibility. This measure consists of four different trials that differentiate between word reading, color naming, inhibitory control and cognitive flexibility (Delis et al., 2001). The participants were asked to name the colors of square blocks on a page (Trial 1, Color Naming), read color words that are printed in black (Trial 2, Word Reading), name the color of the ink that the word is printed in (Trial 3, Inhibition), and to switch between naming the color of the ink that the word is printed in and read the actual word based on a rule (Trial 4, Inhibition/Switching). Each trial was preceded by a sample trial to ensure that the examinee understood instructions. Trial 4 measures cognitive flexibility, as the examinee is required to switch between inhibitory and non-inhibitory responses. The amount of time that it took to complete the task was transformed into z-scores based on normative data. Notably, two other cognitive flexibility measures (i.e., Number-Letter Sequencing Trial of the Trail Making Test, Category Switching trial of the Verbal Fluency Test) were given as part of the comprehensive neuropsychological battery. The Inhibition/Switching trial from the Color Word Inference Test was chosen to represent cognitive flexibility in this study based on its superior psychometric properties out of all three measures. The test-retest correlation for the Inhibition/Switching trial of the Color-Word Interference Test is in the moderate range (r = 0.65), whereas the test-retest correlation for the other two subtests are substantially weaker (r = 0.38 for the Number-Letter Sequencing trial of the Trail Making Test, r = 0.36 for the Category Switching trial of the Verbal Fluency test). In addition, the internal consistency of the Inhibition/Switching trial is moderate to high across age groups (range = 0.52 to 0.8) (Delis et al., 2001).

Measure of cumulative treatment and neurological risk factors

The Neurological Predictor Scale (NPS) incorporates different tumor treatments and neurological risk factors into one score. This measure is used to examine how cumulative risk factors affect outcomes particularly in heterogeneous samples in which survivors vary with respect to treatments and other neurological risk factors (Micklewright et al., 2008). The NPS incorporates information about radiotherapy, chemotherapy, neurosurgery treatments, and the presence of other common neurological risk factors such as endocrine dysfunction, hydrocephalus, and seizure medications into a score that ranges from 0 (lowest level of risk) to 11 (highest level of risk) (Micklewright et al., 2008). Studies have documented the reliability and validity in childhood survivors (McCurdy et al., 2016a; Micklewright et al., 2008; Papazoglou et al., 2008). This measure is significantly associated with intelligence, adaptive functioning, processing speed, working memory and attention over and above each individual risk factor (King and Na, 2016; Taiwo et al., 2017).

Neuroimaging measurement

Neuroimaging data was acquired using a 3 T Siemens Trio MRI scanner using a 12-channel head coil. Participants were outfitted with protective earplugs to reduce scanner noise. Diffusion tensor imaging data was acquired using a single-shot spin echo echo-planar imaging (EPI) sequence with 30 gradient directions and the following acquisition parameters: repetition time (TR) = 7700 ms; echo time (TE) = 90 ms; b = 1000 s/mm2; acquisition matrix = 204 × 204; voxel size = 2.0 × 2.0 × 2.0 mm, 60 contiguous axial slices and scan time = 8 min 22 s. High-resolution T1-weighted structural images were also acquired by collecting 176 contiguous sagittal slices using a three-dimensional magnetization prepared rapid gradient echo imaging (3D MPRAGE) sequence with the following parameters: repetition time (TR) = 2250 ms; inversion time (TI) = 850 ms; echo time (TE) = 3.98 ms; field of view (FOV) = 256 mm; acquisition matrix = 256 × 256; voxel size = 1.0 × 1.0 × 1.0 mm; slice thickness = 1.0 mm; flip angle = 9°. A field map was also recorded with a gradient echo sequence with the parameters of repetition time (TR) = 488 ms; echo time 1 (TE 1) = 4.92 ms; echo time 2 (TE 2) = 7.38 ms; voxel size = 3.0 × 3.0 × 3.0 mm; field of view (FOV) = 204 mm; slice thickness = 3.0 mm; 40 slices; flip angle = 60° to measure field inhomogeneities and compensate for geometrical distortions that result from standard EPI sequences.

Procedure for image analysis

The pipeline for the study included a series of steps for preprocessing, tractography and network construction, detailed below and in Fig. 1.

Fig. 1

Data processing pipeline. A. Diffusion tensors were calculated, FA maps were generated and deterministic tractography was conducted. B. Diffusion and T1 images were co-registered, and the co-registered image was registered to standard space. C. These transformation matrices were combined, inversed and applied to the AAL to yield a parcellation in native diffusion space for each participant. D. A weighted adjacency matrix was created where edges were defined as the average FA of all the voxels along streamlines linking two nodes. E. The Brain Connectivity Toolbox was used to calculate topological properties.

Preprocessing

Diffusion-weighted images underwent visual inspection for distortion, artifact, or clear movement that may render the image unusable for analysis. The images that passed inspection underwent correction for eddy current distortion and subject movement using the “eddy” tool from FSL (Andersson and Sotiropoulos, 2016) and were skull-stripped using the Brain Extraction Tool (Smith, 2002). Results of skull-stripping were manually assessed for quality for all participants.

The estimated translational and rotational displacement for each frame (compared to the frame that immediately preceded it) was quantified in the x, y, and z axes and summarized into one motion metric for each individual as outlined in Power et al. (2012) with the following empirical formula: FD = |∆d| + |∆d| + |∆d| + |∆α| + |∆β| + |∆γ|, where ∆d = d( − d (i.e., the level of translational displacement from one frame to the previous frame in the x-axis) and so on for each of the other parameters. Given that differing levels of motion between groups can have a systematic impact on results, independent samples t-tests and correlations were conducted to determine whether motion represented a confound in the analyses. It was determined a priori to use motion as a covariate for analyses that compared between the two groups if the level of motion differed significantly between groups and correlated significantly with the outcome measure.

The “epi_reg” tool was used to co-register diffusion images to T1-weighted images while correcting for EPI distortions using the fieldmap using the “fsl_prepare_fieldmap” tool in FSL (Jenkinson et al., 2002; Jenkinson and Smith, 2001). The co-registered image was registered to a high resolution standard space using the “auto_warp” command in AFNI which pairs an affine transform (with 12 degrees of freedom) and a volume-based nonlinear transform. The nonlinear registration involved dividing the source image into shrinking and overlapping 3D patches and warping the source image to the template image in incremental steps based on Hermite cubic basic functions. Results of the registration were manually assessed for quality for all participants. Finally, diffusion tensors were calculated and FA maps were generated using FSL's “dtifit” tool.

Tractography and network construction

Deterministic tractography was performed using the Diffusion Toolkit in PANDA, a MATLAB toolbox for pipeline processing of diffusion MRI images (Wang et al., 2007). Whole brain tractography was conducted by placing a seed in all white matter voxels and linearly propagating lines from each seed based on the principal direction of the tensor in that voxel. Each line was propagated by 0.25 mm to the next point in space, at which point the process was repeated. Each of these streamlines was terminated when FA < 0.15 or when the angular deviation from paths was >55° to prevent streamlines from looping back. All possible streamlines were constructed from each seed region. The FA threshold of 0.15 was used as one of the termination criterion as prior research has shown that survivors have overall lower white matter integrity when compared to age-matched controls. The cingulum in the cingulate gyrus part was visualized for several participants using Trackvis based on ROI protocols from prior research (Wakana et al., 2007) to ensure that the whole brain tractography followed the trajectory of a long distance white matter tract (Supplementary).

The Automated Anatomical Labeling Atlas (AAL) was used as the parcellation scheme to indicate nodes of interest for this study (Tzourio-Mazoyer et al., 2002). The AAL atlas divides the brain into 47 cortical volumes of interest in each hemisphere and 26 subcortical regions. Each region in the atlas defined a node in the network analyses, while the average FA between each node pair represented edges in the network analysis. The AAL atlas was transformed to yield a parcellation scheme in native diffusion space; parcellation results were manually assessed for quality for each participant.

The transformed AAL atlas was used to filter the whole brain file to only include streamlines that passed through each node pair. Streamlines with two end-points within the masks of each given node pair were considered to link the two nodes. The average FA of all the voxels along streamlines linking two nodes were defined as the edge weight value for that node pair. An adjacency matrix was constructed where each node was represented in rows and columns and edge values were entered into cells of the intersecting row and column of the corresponding node pair. Finally, the connectivity matrix was averaged across the entire sample and used to generate a binarized mask that included all cells with non-zero values. This mask was then applied to each individual (regardless of group membership), to control for the differences due to the differences in connectivity patterns.

Network properties

The Brain Connectivity Toolbox was used to calculate the topological properties of each network. Measures of properties of nodes and properties representing network integration and segregation were calculated. A short description of the relevant metrics is provided below; more detailed mathematical definitions can be found in Rubinov and Sporns (2010).

Density: Density is a basic characteristic of the network and describes how many existing edges there are in the network out of the number of total possible edges. Methodological studies have demonstrated that other network metrics change as a result of density rather than the topological properties of the network (van Straaten and Stam, 2013). As such, proportional thresholding procedures were used to account for differences in network densities before further analyses. This procedure preserves the same proportion of the strongest edge weights across all individuals.

Global efficiency: Global efficiency is a measure of integration that reflects a characteristic of the overall network. It is calculated as the inverse of the path length, which is defined as the average of the fewest number of edges between all node pairs in the network. A network with a high global efficiency suggests high capacity for parallel processing and thus higher levels of global processing.

Clustering coefficient: The clustering coefficient is a measure of segregation that represents the probability that the neighbors of a node are also connected to each other in the form of a triangle. A node with high clustering coefficient suggests high levels of local processing in that node. The clustering coefficient across all nodes is averaged for an overall measure of segregation in the structural network.

Modularity: Modularity is another measure of segregation and is defined as the existence of communities that have more connections with one another (i.e., high number of within-group links) than is expected in a random model. High modularity values suggest the existence of communities of nodes that have specialized functions.

Betweenness Centrality (BC): BC is a measure of centrality that is calculated as the number of shortest paths that must pass through that node. A node with a high BC suggests that the node is important in the overall network and has a large influence on the transfer of information throughout the overall network.

Hub Disruption Index (HDI): The hub disruption index was calculated from the BC values of all nodes in the network (Termenon et al., 2016). The HDI is defined as the slope of the best-fit line through a plot where the x-axis represents the average BC for each node in the healthy control group, and the y axis represents the difference between the BC for each node between the survivor and the average healthy group. A high negative slope that passes through the x-axis suggests significant and preferential damage to the hubs.

Statistical analyses

Adjustments were made for multiple comparisons to reduce the potential for Type I error. Since there were a priori hypotheses and planned comparisons derived from prior existing literature for each aim, each aim was considered a unique and independent question. Accordingly, we limited our adjustment for multiple comparisons to the number of graph theory metrics that were tested for each aim. As there were four metrics analyzed for each hypothesis, results were considered significant at a p-level equal to or below 0.0125 (i.e., p ≤ .05/4). For the first aim, independent two-sample t-tests were conducted on global efficiency, average clustering coefficient, and modularity to test whether the two groups differed on these metrics. A one-sample t-test was conducted on the average hub disruption values among the survivors to test whether the index was significantly lower than zero. Bivariate Pearson correlations were conducted to test whether younger age at diagnosis, longer time since diagnosis and cumulative neurological risk were associated with more disruptions to white matter network topology. Pearson correlations were used as all variables used for correlation analyses were collected on a continuous or interval scale, and passed assumptions (i.e., no outliers, approximately normal distribution).

To test whether properties of the network related to cognitive flexibility for the second aim, bivariate Pearson correlations were conducted on metrics from the network with age-normed z-scores from the Inhibition/Switching Trial of the Color-Word Interference Test. The second hypothesis of aim 2 was that the differences in white matter network topology would mediate the cognitive differences between survivor and control groups. Given that prior literature has demonstrated that the global efficiency of structural brain networks is most robustly related to executive functions in several different clinical populations including patients with traumatic brain injury, adolescents with congenital heart disease, and adults with mild cognitive impairment, global efficiency was proposed to represent the mediator variable in the mediation analyses (Berlot et al., 2016; Caeyenberghs et al., 2014; Panigrahy et al., 2015).

The SPSS “indirect” script was used to test the mediation model with group membership (survivors vs. controls) as the independent variable, the Inhibition/Switching Trial performance as the dependent variable, and the graph theory metric as the hypothesized mediator. Given the relatively small sample size and the concerns of the Baron and Kenny (1986) model and Sobel test for detecting effect sizes in small samples, bootstrapping was employed with 10,000 samples (Preacher and Hayes, 2004). An effect was deemed significant if the resulting 95% confidence interval of the indirect effect of the independent variable on the dependent variable did not include zero. As this approach can increase the likelihood of Type I error, a test of joint significance was also conducted; if the paths of the regression between the independent variable and the hypothesized mediator (path ‘a’), as well as the regression between the hypothesized mediator and the dependent variable (path ‘b’) were significant, then the indirect effect was also considered statistically significant.

For the third hypothesis, the SPSS “indirect” script and test of joint significance methods were used to test whether the relationship between the NPS score and cognitive flexibility would be mediated by the global efficiency in survivors.

Results

Motion

The mean average displacement for each frame was 0.60 mm (SD = 0.16) for controls and 0.65 mm (SD = 0.16) for survivors. Average displacement did not differ significantly between the two groups, t(74) = 1.25, p = .22, d = 0.29. Motion did not correlate significantly with any of the graph theory metrics or cognitive flexibility performance (p > .05, see Supplemental). Given that motion did not vary between the two groups and did not relate to the dependent variables in the study, motion was not considered a confound and was not used as a covariate for analyses.

Density

The means and standard deviations of the raw, unmasked and mean connectivity matrices are presented in Fig. 2.

Fig. 2

A. Mean and standard deviations of the raw weighted connectivity matrices in control (left) and survivor (right) groups. B. Binarized matrix based on the combined matrices across the entire sample (non-zero mean FAs are indicated in white) C. Scatterplot of mean FA in controls and mean FA in survivors for each edge. The black line indicates a perfect linear correlation between mean FA in controls and survivors for each edge; dots above the black line indicate connections in which controls have higher FAs than survivors. Controls have higher FAs than survivors in 71% of the edges at connections of medium to high FAs, defined as FA > 0.3 (arrows in the center and upper right of the figure), whereas controls have higher FAs than survivors in 56% of the edges at connections of low FA.

Two-sample t-tests were conducted to determine whether network density differed between the two groups. Average density in the healthy controls (M = 0.29, SD = 0.03) was significantly higher than the average density in survivors (M = 0.27, SD = 0.03), t(74) = 2.44, p = .02, d = 0.7. Given that differences in density can drive differences in graph theory metrics that may not reflect differences in topology, proportional thresholding was used based on the average density across the entire sample (density = 0.28), where the strongest 28% of the connections in each network were preserved for both groups. To test whether proportional thresholding may have unduly affected results, the same analyses were run on networks without any thresholding and with stricter proportional thresholding (density = 0.20); see Supplemental for results from different thresholding schemes. As the results did not appreciably change, results presented are based on proportional thresholding using the average density across the sample.

Differences in graph theory metrics between groups

Consistent with hypotheses, global efficiency and average clustering coefficient were higher in controls compared to survivors (see Table 2). Modularity did not differ significantly between the two groups. The HDI was significantly different from zero, which was also consistent with hypotheses, and indicated preferential damage to hub regions. A graphical representation of the mean BC values for controls and survivors, as well as the graphs used to derive the hub disruption index are represented in Fig. 3. Survivors showed lower BC values in regions that are hubs for healthy controls, indicating compromise of regions that are of high importance to brain networks. Further, the negative slope in Fig. 3B confirms that there are larger differences in measures of centrality in hub regions as compared to other nodes.

Table 2

Graph theory metrics in survivors and healthy controls.

Measure	Controls (n = 38)		Survivors (n = 38)		df	t	p	Cohen's d
	M	SD	M	SD
Global efficiency	0.31	0.014	0.29	0.019	74	3.67	0.000	1.20
Avg. clustering coefficient	0.27	0.013	0.26	0.015	74	2.82	0.006	0.71
Modularity	0.25	0.04	0.26	0.05	74	−0.30	0.762	0.22
Hub disruption index			−0.07	0.14	37	−3.18	0.003	0.50

Fig. 3

A. Betweenness centrality (BC) bar plots derived from the control (top) and survivors (bottom) groups. The order for both maps was based on the BC values in controls. The solid vertical bars and error bars represent group BC means and standard deviations, respectively. The solid horizontal lines represent the mean BC across all brain regions in each group, while the dashed horizontal lines correspond to the mean BC plus one standard deviation. The light grey bars in both maps represent regions where BC is higher than the mean plus one standard deviation across all brain regions. B. Graph of the best fit lines for each survivor, where the x-axis represents the mean BC of nodes in the control group, and the y-axis represents the subtraction of the mean BC in the control group from the BC in each survivor for each node. The hub disruption index for each survivor is equal to the slope of the best fit line; a negative slope indicates preferential damage to hub regions.

Correlations with risk factors and cognitive flexibility

Higher scores on the NPS were associated with lower global efficiency and lower average clustering coefficient (Table 3, Fig. 4) after correcting for multiple comparisons.

Table 3

Correlations between risk factors and graph theory metrics (n = 38).

Measure	Graph theory metric
	Global efficiency	Avg. clustering coefficient	Modularity	Hub disruption index
Age of survivor at diagnosis	−0.029	0.06	0.11	0.12
Time between diagnosis and exam	−0.22	−0.27	0.008	−0.13
Neurological Predictor Scale	−0.61*	−0.65*	0.23	−0.08

Note. *p < .0125 (significant after Bonferroni corrections for multiple comparisons).

Fig. 4

Scatterplots of correlations between NPS Score and graph theory metrics in survivors.

After Bonferroni corrections for multiple comparisons, global efficiency was significantly correlated with cognitive flexibility, as measured by performance on the Inhibition/Switching trial of the Color Word Interference Test, r(74) = 0.40, p < .05. Average clustering coefficient was also significantly associated with cognitive flexibility, r(74) = 0.35, p < .05 (corrected). Modularity and HDI were not statistically correlated with cognitive flexibility.

We hypothesized that differences in brain network properties would mediate the difference in cognitive flexibility between survivor and control groups. Given that global efficiency was most correlated with scores on the Inhibition/Switching Trial of the Color-Word Interference Test and had the highest effect sizes for group differences, global efficiency was used as the mediator in the model. The confidence interval for the indirect path (path c’) did not include 0 and both paths a and b in the model were significant, indicating that global efficiency mediated the differences in cognitive flexibility performance between the two groups (Fig. 5). Notably, the direct effect of the independent variable group membership on cognitive flexibility was not statistically significant. Although traditional approaches to mediation analyses require a significant direct effect of the independent variable on the outcome variable to test for mediation, more modern statistical perspectives posit that significant indirect effects through mediators do not depend on the presence of statistically significant direct effects, especially within the context of a theoretically meaningful model (Hayes, 2009).

Fig. 5

Global efficiency mediates cognitive flexibility differences between groups.

Similarly, the confidence interval for the indirect effect of NPS score on cognitive flexibility (path c’) did not include zero, and both paths and b of the model were statistically significant (Fig. 6). These results suggest that the association between cumulative neurological risk and cognitive flexibility was explained by the global efficiency of the structural network.

Fig. 6

Global efficiency mediates the relationship between NPS and cognitive flexibility in survivors.

Discussion

The results of this study indicated that global efficiency and average clustering coefficient of white matter networks were reduced in survivors of pediatric brain tumors compared to healthy peers matched by age, gender, handedness, and socioeconomic status. There was also evidence for preferential impact to hub regions. Further, lower global efficiency and lower average clustering coefficient were associated with higher cumulative neurological risk and poorer performance on behavioral measures of cognitive flexibility. Indeed, global efficiency mediated differences in cognitive flexibility performance between survivors and healthy peers, as well as the relationship between cumulative neurological risk and cognitive flexibility performance. These results suggest that structural networks are altered in adult survivors of pediatric brain tumors and that topological features of these networks explain differences in cognitive flexibility performance. These results are highly consistent with findings from studies conducted in other clinical groups such as TBI, stroke, epilepsy, and congenital heart disease, which have shown disruptions in measures of segregation, integration, and centrality when compared to healthy adults. Prior studies also have consistently shown that metrics describing the integrity of the network significantly relate to behavior and the degree of impairment.

Global efficiency, a measure of global integration, is hypothesized to reflect the capacity of networks to allow efficient processing of information from distributed regions of the brain. The clustering coefficient is a measure of segregation that represents high levels of local processing. Brain networks of healthy individuals are associated with a balance of segregation and integration in the brain so that information can be efficiently transferred across and between brain networks while still maintaining low biological costs. Studies of healthy developing brains indicate that brain networks undergo highly dynamic changes from infanthood to late adolescence; these networks change from relatively random configurations to networks that optimize the balance between information segregation and integration. These changes support cognitive and behavioral developments (Baum et al., 2017; Cao et al., 2017; Chen et al., 2013). The rapid changes occurring in brain networks during development also render the brain more vulnerable to neurological insults. Survivors of pediatric brain tumors experience disruptions during these critical timeframes when structural and functional networks are actively being optimized for efficiency. The results of this study suggest that these network alterations persist when survivors have grown into adulthood. These alterations also may represent the neural underpinnings of behavioral outcomes, as reduced global efficiency and reduced average clustering coefficients are associated with poorer cognitive flexibility.

Inconsistent with hypotheses, modularity of the structural networks was not different between groups and was not significantly related to measures of cognitive flexibility. Lack of expected findings may be due to the parcellation scheme. As the AAL's nodes are defined by anatomical boundaries, it is possible that nodes contain subregions that are heterogeneous with regard to function and cytoarchitecture and thus the parcellation may not have the spatial resolution to be sensitive to differences in modularity. Alternatively, recovery may prioritize modularity of structural networks over other types of global and local connectivity. Computational modeling studies that incorporate plasticity into their network recovery models have reported that modularity recovers over time after a lesion in the network (Stam et al., 2010). The cross-sectional design of this study precludes examining how modularity changes over time on an individual level. As such, longitudinal studies will be crucial to establish whether acute and subacute stages of injury are associated with changes in modularity and whether modularity recovers as a function of time. It is also possible that differences in modularity might exist in functional networks during the chronic phase of injury that are not reflected in the structural measures used for this study. Although functional networks are constrained by underlying structural architecture, functional networks also exhibit more flexibility and can adapt quickly to environmental demands by reorganizing, coordinating and mobilizing different regions across the brain (Fischer et al., 2014). Multi-modal network approaches will be helpful in understanding how the brain recovers over time and how these changes contribute to behavior in adult survivors of pediatric brain tumors.

Our results also showed that survivors had lower BC values in regions that are hubs for healthy controls, indicating compromise of regions that are of high importance to brain networks. Further, the HDI was significantly different from zero, further confirming that there were larger differences in measures of centrality in hub regions as compared to other nodes. This is consistent with prior literature indicating preferential effects to hub regions across a variety of brain disorders (Crossley et al., 2014). Indeed, it should be noted that preferential effects to hub regions were found in a sample of survivors who were heterogeneous with regard to tumor type and location. This suggests that there may be a level of specificity to the regions in the brain that are most at risk across long term survivors of varying types and locations of brain tumors.

It is important to interpret the findings of this study within the context of the methods used, as graph theory approaches to study the brain are still in development and the biological significance of these metrics is still under investigation (He and Evans, 2010). The lack of a gold standard for processing or conducting complex network analyses in clinical populations have potential implications for findings, as methodological choices on parcellation schemes, edge definitions and thresholding procedures can significantly impact the network metrics under investigation (van Wijk et al., 2010). Given the lack of a gold standard, this study used similar methods, parcellation schemes, edge definitions, and thresholding procedures as other research studies investigating structural network properties in clinical populations. Notably, research has shown that graphs with different densities are difficult to compare directly because density differences in networks can result in significant differences in graph theory metrics even when networks share the same topological organization. This makes it challenging to compare networks between clinical groups and healthy same-aged peers when density may naturally change as the function of the disorder itself. There are different methods that have been proposed to deal with comparing networks of different densities but there is no one satisfactory way to control for this issue completely, due to the fact that modeling the exact impact of density on different graph theory metrics depends on knowing the underlying topology a priori, which is not possible in empirical studies of clinical populations (van Wijk et al., 2010). For instance, some studies use density as a covariate to test whether graph theory metrics remain different between groups after removing the variance associated with density. This approach, however, does not fully control for the issue when density does not share linear relationships with graph theory metrics (Caeyenberghs et al., 2012). Proportional thresholding, which was the method used in this study, uses a cutoff such that the same percentage of edges are enforced for everyone's networks. This method, however, may lead to modifications of the network by ignoring significant connections in controls or by enforcing weaker connections in the clinical group (Drakesmith et al., 2015). In addition, given that controls had higher densities, using any cutoff automatically affects the control group more than the survivor group. To test whether proportional thresholding may have unduly affected the results, the same analyses were run on networks without any thresholding. The results were very similar (see Supplemental); global efficiency and clustering coefficients remained lower in survivors compared to controls, were significantly correlated with cognitive flexibility and cumulative risk, and mediated the difference in cognitive flexibility performance between the two groups. The same sets of analyses were also conducted on a stricter thresholding level (density = 0.2). Group differences in global efficiency remained, and both the average clustering coefficient and global efficiency correlated significantly with cognitive flexibility and cumulative neurological risk. One notable difference in the analyses with the strictest thresholding procedure was that the clustering coefficient was no longer significantly different due to a larger decrease in clustering coefficient values in controls compared to survivors. This possibly suggests that the “weaker” edges removed due to the strict thresholding procedure contributes more heavily to a clustered and segregated network in healthy individuals. Notably, much of the research on the topological features of brains in healthy and clinical groups have largely assumed that “stronger connections are better” and have concentrated their efforts on understanding the nature of these strongest connections by using thresholding procedures across that preserve the strongest connections and remove the weakest connections when examining networks. However, recent research has recognized the importance of these weaker connections in explaining individual variability and symptom presentation in clinical disorders (Bassett et al., 2012; Santarnecchi et al., 2014). These additional analyses suggest that the results of this study are not fully attributable to densities of networks but does warrant continued investigation of the impact of density on topological features of networks in adult survivors of pediatric brain tumor.

The findings from this study should be considered within the context of the limitations. First, both survivor and control groups were self-selected, which may have skewed the sample to higher functioning survivors who may have had the time and means to transport themselves to the study site or lower functioning survivors with more cognitive concerns. Due to these factors, selection bias may limit the generalizability of the conclusions. In addition, the control group was more ethnically diverse than the survivor group. Post hoc analyses, however, confirmed that group differences in topology remained after adjusting for ethnicity. Furthermore, the correlations between the network characteristics and neurological risk, as well as cognitive flexibility, remained after adjusting for ethnicity. Related to neuroimaging methodology, there are established limitations with the diffusion weighted imaging parameters and deterministic tractography methods used in this study. Research has established that false positive and false negative streamlines can occur due to signal noise, partial volume effects and complex fiber architecture within voxels (Jbabdi et al., 2015) and that these methods can lead to bias in FA estimates especially in long white matter tracts (Oouchi et al., 2007). Future studies will need to employ more advanced diffusion imaging models, such as high angular resolution diffusion imaging or diffusion spectrum magnetic resonance imaging to more accurately track complex fiber architecture in regions where more complex fiber architectures exist.

There are also several strengths that are worthy of note. Like research studies of other clinical populations, the sample studied in this group is heterogeneous with regard to tumor type, tumor location, and level of neurological risk. This heterogeneity allowed for increased variance to explore the relationship between variables of interest and functional outcomes. Nevertheless, to test whether the heterogeneity of the study sample may have obscured or artificially created a relationship between variables, a subgroup analysis was conducted in the subsample of survivors with tumors located in the posterior fossa (n = 26, see Supplemental for results of analyses). The findings from the subsample were highly consistent with the findings from the entire group; global efficiency and average clustering coefficient of white matter networks were reduced in survivors compared to healthy controls, with preferential impact to hub regions. Further, lower global efficiency continued to be associated with higher cumulative neurological risk and poorer performance on behavioral measures of cognitive flexibility. This subgroup analysis suggests that the findings from the study are not strictly due to the heterogeneity in tumor locations and further strengthens the conclusions drawn from the study.

Other strengths of the study include a respectable sample size of an understudied population, as well as the use of an age- and gender-matched control group that allowed for comparisons between the clinical group and healthy same aged peers and to determine whether structural topology was altered in the survivor group and underpinned differences in cognitive performance between groups.

To our knowledge, this is the first study to use graph theory to explore the topological properties of white matter networks in survivors of brain tumors. Understanding structural topology lays the groundwork for exploring functional network organization, as functional networks are shaped and constrained to a certain extent by the underlying structure (Cao et al., 2017). Understanding the flexibility and diversity of functional network organizations within the constraints of anatomical connectivity can provide important insights into the nature of brain repair, recovery, and function following a neurological insult. Metrics derived from structural and functional brain networks have been used as a biomarker in clinical groups such as patients with temporal lobe epilepsy to predict patients who will have better outcomes after surgery (Bonilha et al., 2015; Ji et al., 2015). Studies of this kind suggest that graph theory may even have some utility in clinical settings to guide neurosurgical planning to avoid neurological deficit, predict the efficacy of treatments, and identify patients who are at risk for poor outcomes (Castellanos et al., 2013; Petrella, 2011). As this study was cross sectional and thus could not provide information about causation, or about changes that occur over time, longitudinal and prospective studies will be crucial to establish time frames and causes of change. For instance, to demonstrate that structural network topology truly plays a causative role in explaining cognitive flexibility performance outcomes, it will be important to establish that structural network changes precede behavioral changes. Although longitudinal work is necessary in larger samples to establish that graph theory metrics have clinical value in survivors, the findings from this study suggest the potential clinical relevance of understanding white matter network properties of the brain. The current study demonstrated the efficacy of graph theory when examining cognitive flexibility and white matter networks in adult survivors of childhood brain tumors and highlighted the importance of future longitudinal studies to characterize long-term outcomes in personalized medicine.

Declarations of interest

None.