Abnormal structural connectivity in the brain networks of children with hydrocephalus.

Increased intracranial pressure and ventriculomegaly in children with hydrocephalus are known to have adverse effects on white matter structure. This study seeks to investigate the impact of hydrocephalus on topological features of brain networks in children. The goal was to investigate structural network connectivity, at both global and regional levels, in the brains in children with hydrocephalus using graph theory analysis and diffusion tensor tractography. Three groups of children were included in the study (29 normally developing controls, 9 preoperative hydrocephalus patients, and 17 postoperative hydrocephalus patients). Graph theory analysis was applied to calculate the global network measures including small-worldness, normalized clustering coefficients, normalized characteristic path length, global efficiency, and modularity. Abnormalities in regional network parameters, including nodal degree, local efficiency, clustering coefficient, and betweenness centrality, were also compared between the two patients groups (separately) and the controls using two tailed t-test at significance level of p < 0.05 (corrected for multiple comparison). Children with hydrocephalus in both the preoperative and postoperative groups were found to have significantly lower small-worldness and lower normalized clustering coefficient than controls. Children with hydrocephalus in the postoperative group were also found to have significantly lower normalized characteristic path length and lower modularity. At regional level, significant group differences (or differences at trend level) in regional network measures were found between hydrocephalus patients and the controls in a series of brain regions including the medial occipital gyrus, medial frontal gyrus, thalamus, cingulate gyrus, lingual gyrus, rectal gyrus, caudate, cuneus, and insular. Our data showed that structural connectivity analysis using graph theory and diffusion tensor tractography is sensitive to detect abnormalities of brain network connectivity associated with hydrocephalus at both global and regional levels, thus providing a new avenue for potential diagnosis and prognosis tool for children with hydrocephalus.

Introduction

Normal brain function requires the integrity of neuronal function and connectivity, globally and regionally. In HCP, data from both human research and experimental research have suggested that a wide range of WM networks connecting various functionally important cortical and subcortical regions are primary targets for disruption due to enlarged ventricles and/or increased intracranial pressure (Assaf et al., 2006; Hassan et al., 2008; Yuan et al., 2009; Yuan et al., 2010; Scheel et al., 2012; Ginat et al., 2013; Yuan et al., 2013). Although surgery can significantly reduce the mortality and morbidity, HCP remains as an incurable lifelong disorder (Mataro et al., 2001). The damage to neuroanatomy sustained prior to the surgical treatment may remain or continue to progress after the surgery, leading to long term behavioral and neuropsychological deficits in visuospatial skills, visuomotor skills, and a series of other important neurocognitive domains (Erickson et al., 2001; Mataro et al., 2001; Frank et al., 2003; Bakar et al., 2009).

DTI is an advanced neuroimaging technique that can measure in vivo WM structural integrity (Basser and Jones, 2002; Beaulieu, 2002). A growing body of literature has shown the success of DTI in investigating WM structural abnormality in pediatric HCP (Assaf et al., 2006; Hassan et al., 2008; Yuan et al., 2009, 2013; Air et al., 2010; Jang and Kim, 2011; Sheel et al., 2012; Jiang et al., 2013; Rajigopal et al., 2013). A series of regions, e.g., corpus callosum and posterior internal capsule, have been reported to show different degrees and patterns of abnormality in the directionality as well as the magnitude of water diffusion parameters. In these studies, in order to extract DTI parameters, a predefined list of ROIs needs to be delineated based on the hypotheses of the study. However, this approach does not provide a global quantification for the integrity of the entire brain network, nor does it allow for detecting any abnormalities outside the hypothesized brain regions. Hydrocephalus, a neurological disorder with heterogeneous etiologies, often presents with regionally specific damage, and this damage is also expected to extend to wider areas throughout the brain. ROI-based DTI analysis alone is clearly not sufficient to meet all the challenges in studying this patient population.

In recent years, graph theory analysis has emerged as a promising tool that provides information on brain connectivity, structural and functional, at both global and regional levels (Rubinov and Sporns, 2010). The brain is modeled as a network composed of a number of nodes and edges connecting these nodes. The nodes and edges respectively represent individual gray matter regions responsible for various brain functions and white matter fibers responsible for transferring information among these regions. This method of graph theory analysis has been applied in characterizing the developmental trajectory of network connectivity as well as for investigating the disruption of network connectivity in various neurological disorders (Liu et al., 2008; Bernhardt et al., 2011; Kim et al., 2011; Shi et al., 2012; Ottet, 2013). So far, however, no study has used graph theory approach in the analysis of white matter integrity in hydrocephalus brain networks.

In the present study, graph theory is applied to analyze the structural connectivity based on DTI tractography in children with hydrocephalus. The hypothesis of the present study is that the brain network integrity is affected in hydrocephalus patients as reflected in the aberrant topological features at both the global and regional levels. More specifically, we aim to investigate whether brain networks in children with hydrocephalus exhibit small-world properties and to evaluate whether small-worldness, a parameter derived from the graph theoretical analysis, is a sensitive measure for detecting global network alteration in this patient population. Specifically, we hypothesize that the small-worldness, and other global network measures, including normalized clustering coefficient, normalized characteristic path length, global efficiency, and modularity, are abnormal in both hydrocephalus patient groups. In addition, we aim to assess abnormalities of regional network measures, including nodal degree, local efficiency, clustering coefficient, and betweenness centrality, in hydrocephalus patients at both pre-surgery and post-surgery in comparison to normal controls.

Materials and methods

Patients

This was a retrospective analysis with all the data selected from an ongoing prospective neuroimaging project of children with hydrocephalus before and within 1 year after CSF diversionary surgery. At the time of data analysis, 61 normal children and 58 children with hydrocephalus were recruited into the study. It was decided that we would need to (1) exclude datasets from children younger than 11 months because the image normalization did not generate consistent results due to the poor image contrast in very young children; and (2) exclude datasets that presented image artifact because of the programmable valves in the shunts. In addition, five participants (3 controls and two HCP patients) were initially eligible for the study but were also needed to be excluded due to excessive head motion (n = 4) or image distortion (n = 1). Combining these factors, three groups of children were included in the study: Group 1 were controls, n = 29, age range: 13.1–197.8 months, median age 48.7 months; Group 2 were preoperative (PreOp) HCP patients, n = 9, age range: 11.0–194.5 months, median age 38.5 months; and Group 3 were postoperative (PostOp) HCP patients, n = 17, age range: 12.0–207.2 months, median age 41.4 months. No statistically significant difference was found in age between the control group and either one of the two patient groups (p > 0.1). The demographic information for patients is included in Table 1.

Table 1

Study population demographics.

Subject	Gender	Gestational age (wks)	Birth weight (g)	Etiology of HCP	Additional MRI findings	Additional pathologies; neurological/psychological disorder
Subj_01	M	36	4040	Congenital HCP; AS	Tectal dysplasia;cerebellar tonsillar herniation; periventricular leukomalacia; Chiari malformation	Developmentally delayed
Subj_02	M	34	1400	Congenital HCP	None	Global developmental delay;spastic quadriplegic cerebral palsy
Subj_03	M	40	3200	Congenital HCP; AS	CC thinning, WM injury	Developmentally delayed
Subj_04	M	35	Unknown	Communicating HCP	None	None
Subj_05	F	40	2500	Congenital HCP	None	Developmentally delayed
Subj_06	F	34	2000	Congenital HCP	None	None
Subj_07	F	40	2800	Congenital HCP	None	None
Subj_08	M	41	3300	AS	None	Developmentally delayed
Subj_09	M	40	Unknown	IVH; AS	Hemorrhagic venous infarction; cerebral palsy	Developmentally delayed
Subj_10	M	37	2690	Posterior fossa arachnoid cyst	None	Gross motor delays
Subj_11	M	40	2890	Communicating HCP	Chiari 1 (developed after shunt surgery; surgically treated 45 mon after initial shunt surgery)	Transient tic disorder;expressive language disorder
Subj_12	M	39	4000	Obstructive HCP-third ventricular arachnoid cyst	Cystic lesion	Static motor deficits
Subj_13	M	40	Unknown	Communicating HCP	None	Developmentally delayed;language/speech delay
Subj_14	F	38	3800	Congenital HCP	None	None
Subj_15	M	35	2300	Obstructive HCP-tectal lesion	Cystic lesion;cerebral aqueduct lesion	None
Subj_16	F	36	3500	AS	None	None
Subj_17	M	40	3600	AS	None	None
Subj_18	F	Unknown	Unknown	Tectal plate glioma	None	Vitamin D deficiency
Subj_19	M	40	3200	AS	None	None

All the participants were recruited from two hospitals, Cincinnati Children's Hospital (CCHMC) and St. Louis Children's Hospital (SLCH). The study was approved by the Institutional Review Board at both CCHMC and SLCH (Washington University). Families of participants gave written informed consent when enrolled into the study, and children older than 11 years of age provided written assent.

MRI/DTI data acquisition

DTI data were acquired on 1.5 Tesla scanners with a single-shot echo planner imaging sequence at either CCHMC (GE, Signa, GE Healthcare, Milwaukee, Wisconsin) or SLCH (Siemens Avanto, Erlangen, Germany). The sequence specifications were: TR/TE = 9400/93.2 ms; in-plan resolution = 2.5 × 2.5 mm; slice thickness = 2.5 mm; 15 non-colinear diffusion-weighted directions (b = 1000 s/mm2); 1 vol of images with no diffusion sensitization; ASSET or IPAT factor = 2; number of average = 2. Site compatibility was established before the study started and the quality assurance procedures involving both MR phantom and traveling human subject “phantom” were followed rigorously as reported elsewhere (Yuan et al., 2011; Rajagopal et al., 2013; Yuan et al., 2013).

DTI data preprocessing and brain parcellation

The DTI data were corrected for head motion and eddy current artifact using Automated Image Registration, an affine transformation method as described by Woods et al. (1998a,b). The B-matrix was rotated when correcting for subject motion (Leeman and Jones, 2009). To address the recent concerns in the literature about the potential confounding effect of head motion on fiber tracking (Yendiki et al., 2013) and graph theory based connectivity analysis (Power et al., 2012; Satterthwaite et al., 2012; Van Dijk et al., 2012), any datasets that had translational motion that exceeded 1.5 mm and/or rotational motion that exceeded 1 degree were discarded from further analysis. The translation motion in a scan was defined as the median of the frame-by-frame translation values determined by the translation motion in the x, y and z directions. The rotation in a scan was defined as the median of the frame-by-frame rotation values determined as the average of three Euler angles. As described earlier, data of four participants (3 controls, 1 HCP patient) who were initially eligible were discarded because the head motion exceeded 1.5 mm of translation and/or 1 degree of rotation (average of the three Euler angles). In the three groups of participants that were included in the final analysis, no statistically significant difference was found in head motion between the controls and either of the two HCP patient groups.

The DTI metrics maps were calculated with standard technique (Basser and Pierpaoli, 1998). Large deformation diffeomorphic metric mapping, a non-linear transformation (Miller, 2005), was used to normalize the images to MNI space to register with the JHU-DTI-WMPM II atlas (Oishi et al., 2008; Oishi et al., 2009; Djamanakova et al., 2013). DTI maps with multiple contrasts (FA and b0) were used to provide complementary contrast in the normalization. Mask of ventricles in b0 maps was used for HCP patients to improve registration. All the results of the normalization were visually reviewed (W.Y) to confirm accuracy of the procedure (see example in Fig. 1). The parcellation of gray and white matter from the atlas was then inversely transformed back to the subject space and used for parcellating brain into 130 brain regions. From these regions, the 62 cortical and subcortical regions (31 for each hemisphere) were retained for the later fiber counting procedure (Table 2).

Fig. 1

Examples of normalization and parcellation to JHU-DTI-WMPM atlas. (A) Control; (B) HCP patient.

Table 2

Cortical and subcortical brain regions defined in the JHU atlas. 31 regions for each hemisphere.

Region name	Abbreviation	Region name	Abbreviation
Superior parietal gyrus	SPG	Entorhinal area	ENT
Cingulate gyrus	CingG	Superior temporal gyrus	STG
Superior frontal gyrus	SFG	Inferior temporal gyrus	ITG
Middle frontal gyrus	MFG	Middle temporal gyrus	MTG
Inferior frontal gyrus	IFG	Lateral frontoorbital gyrus	LFOG
Precentral gyrus	PrCG	Middle frontoorbital gyrus	MFOG
Postcentral gyrus	PoCG	Supramarginal gyrus	SMG
Angular gyrus	AG	Gyrus rectus	RG
Precuneus	PrCu	Insular	Ins
Cuneus	Cu	Amygdala	Amyg
Lingual gyrus	LG	Hippocampus	Hippo
Fusiform gyrus	Fu	Caudate nucleus	Caud
Parahippocampal gyrus	PHG	Putamen	Put
Superior occipital gyrus	SOG	Thalamus	Thal
Inferior occipital gyrus	IOG	Globus pallidus	GP
Middle occipital gyrus	MOG

In subjects' native space, whole brain fiber tracking was conducted using Diffusion Toolkit/TrackVis (Wang et al., 2007; ) with the entire brain selected as the seed ROI. Tractography was performed based on a deterministic tracking algorithm (FACT: fiber assignment by continuous tracking). An angular threshold of 70° and a FA threshold of 0.15 were applied for the tractography.

Construction of connectivity matrix based on anatomical network

The UCLA Multimodal Connectivity Package (Basset et al., 2011) was used to construct connectivity matrices for graph analysis. As described above, for each subject, brain regions initially parcellated in MNI (Montreal Neurological Institute) space were inversely warped back to subjects' native space. The 62 gray matter regions were included in connectivity analysis and the number of streamlines between each pair of the regions was counted based on the results generated from the whole brain WM fiber tracking. A connectivity matrix (the graph G), was thus constructed, which was a 62 × 62 square matrix with the value in each entry equal to the number of streamlines connecting the corresponding pair of brain regions (Matrix x- and y-indices). The connectivity matrix represents a network that has 62 nodes and undirected edges between the nodes. Fig. 2 shows the schematic diagram for the major steps used to generate the connectivity matrix.

Fig. 2

Flowchart for constructing a structural connectivity matrix. For each subject, DTI data were first preprocessed to minimize effect of head motion and eddy current artifact and then were used to construction of tensors. Whole brain tractography was then performed to generate WM streamlines. The b0 map and FA map (and ventricle masks for HCP patients) were used for registration to the JHU-DTI-WMPM II atlas with the large deformation diffeomorphic metric mapping algorithm. The inverse transformation matrix was used to determine the 62 brain regions in the subject's native space. The number of WM streamlines was calculated for each pair of 62 brain regions. The brain network was then built using the 62 brain regions as nodes and the number of streamlines as edges. The final step was binarizing the initial matrix into the final connectivity matrix with certain threshold. In the present study, the threshold is network wire cost (a.k.a. density) of 0.21.

The original connectivity matrix was binarized based on whether brain regions are connected with WM fibers. If the number of WM streamlines was equal to or greater than 1, the two regions were considered connected and the corresponding entry in the matrix was assigned a value of 1. If no WM streamline was generated in the fiber tracking between two regions, these two regions were considered as not connected and the corresponding entry was assigned a value of 0. In order to ensure that the network topological features were mathematically comparable, a common practice is to compare network features across subjects at a fixed network density level. In this study, the final binarized connectivity matrix was generated with a network density level of 0.21 for all the subjects before the matrix was used for graph analysis. This was the cost level at which all the 62 nodes in the network remained fully connected for all the subjects. Some networks began to have isolated nodes when the density level was set to be under 0.21.

Graph theoretical analysis

Graph theoretical analysis (or Graph Analysis) was performed using the Brain Connectivity Toolbox (Rubinov and Sporns, 2010). Additional in-house Matlab scripts were used to calculate different variables for network connectivity measures. Standard methods were used in the present study, and all the formulas for calculating the network measures can be found in the literature (Rubinov and Sporns, 2010). Only general descriptions are provided here.

We first analyzed small world network properties proposed originally by Watts and Strogatz (1998) using characteristic path length (Lp), network clustering coefficient (C), and small-worldness (σ). Two additional global network properties: global efficiency (Eglob) and modularity (MOD), were also examined.

A structural network G is defined by a set of nodes N and the connections (edges) that link these nodes. For each node, the number of connections linking the node with others is defined as the nodal degree. The distance between two different nodes is defined by the number of edges along the path, the smallest of which is called the shortest path length between these two nodes. The average of all the shortest path length is named the characteristic path length (Lp), a global property that quantifies the level of integration in the network. The network clustering coefficient (C) is the average of nodal clustering coefficient, which is calculated as the ratio between the existing number of edges among the neighbors of a node and the maximum possible number of connections among these neighbors. Contrary to the characteristic path length, network clustering coefficient is considered as a measure to quantify the level of segregation in the entire system. Small worldness (σ) is calculated as the ratio between the normalized clustering coefficient (γ) and the normalized characteristic path length (λ). In order to normalize the two variables, a null random network is constructed in such a way so that both the degree and degree distribution are preserved as compared to the real brain network under investigation. In the present study, the rewiring was repeated 1000 times and the average value of the null random network (Crand and Lprand) was used as the bases to normalize C and Lp (γ = C/Crand; λ = Lp/Lprand). The small-worldness was therefore calculated as the ratio between these two normalized variable (σ = γ/λ = (C/Crand)/(Lp/Lprand)). In theory, the small-worldness variable quantifies the balance between the segregation and integration for the information processing and communication in the system. Typically a network is considered to have small world feature if it satisfies the following conditions: γ ≫ 1, λ ≈ 1, and σ ≫ 1 (Watts and Strogatz, 1998; Achard et al., 2006; Humphries et al., 2006).

The second part of the analysis was centered on the regional network measures in the patient network in comparison to that in the normal controls. Group difference was tested for each node in the network in four nodal network properties: degree, local efficiency, betweenness centrality, and clustering coefficient. The nodal degree is the most basic and fundamental variable in network study. It is the number of connections linking the node with others. For a certain node, the fraction of shortest path between all other pairs of nodes in the network that actually pass through the node is called nodal betweenness centrality, which is a variable that reflects the importance of the node in information transfer with other nodes. Local efficiency for a certain node is the inverse of the average of the shortest path length in the subgraph defined as the set of nodes that are the neighbors of the node of interest. It is often used to reflect system efficiency, redundancy and tolerance to attack.

Statistical analysis

Our initial testing showed that some network measures (global and regional) had age dependence in children in the control group, which was in line with the developmental changes in children reported in the literature (Hagmann et al., 2010; Dennis et al., 2013). To account for the age effect, a linear regression analysis was applied for all global and regional measures and the residual values (the difference between the data and the fitted curve) were used for all the comparisons in the present study. Group differences in all network connectivity variables were tested with two tailed t-test at significance level of p < 0.05 with multiple comparison corrected using false discovery rate (FDR) method. At global level, the FDR method was applied for each patient group separately to correct for potential false positive findings in the five global network measures. At regional level, the FDR correction was also applied for each patient group separately but across all four regional network measures and the 62 network nodes. In addition to the FDR method, we also performed permutation test (Good, 2005) to compare the PreOp group to the Control group in order to minimize potential spurious findings due to the large difference in sample size in these two groups (n = 9 vs. n = 29). This is a method that gives a simple way to compute the sampling distribution for any test statistics, under the null hypothesis that there is no difference between groups in the measured outcome. To conduct this test, we randomly shuffled the groups and generated one million permutations (a random sample from all possible permutation). The ranking of the t statistic from the original data gives a p-value for the group difference. The findings based on the permutation approach (see Supplemental Tables S1 and S2) are consistent with the findings based on the FDR method, mostly with lower p-value from the permutation method. Therefore, the overall results will still be based on the FDR approach to keep the statistical approach consistent throughout the study and to remain conservative in reporting our findings considering the small sample size in the patient group.

Six children with HCP in the PreOp group had both pre- and postop imaging and thus were also included in the PostOp group. Initial longitudinal analysis using paired t-test did not show any significant pre- vs. post-op change, and thus the results were not included in the final analysis. All the group differences presented in the present study were tested cross-sectionally between the control group and the two patient groups separately. Although some of the data for the two patient groups were summarized together because the results were largely similar (in comparison to controls), no intention was made to infer or imply any longitudinal conclusion based on the data presented in this study.

Results

Abnormal global network measures in children with HCP

As shown in Fig. 3, children with HCP in the PreOp group were found to have significantly lower small-worldness (p = 0.005) and normalized clustering coefficient (p = 0.030) when compared to the controls (Table 3). No statistically significant difference was found in normalized characteristic path length, global efficiency, or modularity in children with HCP in the PreOp group when compared to the control group. It was found that all the significant group differences based on the corrected p-values using FDR approach remained statistically significant based on the permutation approach, as shown in Table S1.

Fig. 3

Comparison of global network topology. (A) Small-worldness; (B) normalized clustering coefficient; (C) modularity.

Table 3

Global network measures (all network values are residuals based on linear regression to account for age factor; all p values are FDR corrected).

Global network measures	CTLn = 29	Preop HCPn = 9				Postop HCPn = 17
(Residual value)	Mean ± std	Mean ± std	df	t	p	Mean ± std	df	t	P
γ	0.0008 ± 0.1159	−0.1539 ± 0.1139	36	−3.51	0.030*	−0.1975 ± 0.2073	44	−4.17	0.001*
λ	0.0027 ± 0.0230	−0.0107 ± 0.0192	36	−1.58	ns	−0.0194 ± 0.0362	44	−2.56	0.018*
σ	−0.0021 ± 0.0774	−0.1200 ± 0.1090	36	−3.62	0.005*	−0.1498 ± 0.1379	44	−4.67	0.001*
Eglob	0.0000 ± 0.0059	0.0042 ± 0.0059	36	1.90	ns	0.0012 ± 0.0240	44	0.33	ns
MOD	0.0003 ± 0.0208	−0.0107 ± 0.0191	36	−1.41	ns	−0.0206 ± 0.0248	44	−3.05	0.006*

Note: γ = normalized clustering coefficient; λ = normalized characteristic path length; σ = small-worldness; Eglob = global efficiency; MOD = modularity. ns = not significant.

Children with HCP in the PostOp group were also found to have significantly lower small-worldness (p = 0.001) and lower normalized clustering coefficients (p = 0.001). In addition, children with HCP in the PostOp group were found to have significantly lower normalized characteristic path length (p = 0.018) and lower modularity (p = 0.007) when compared to the control group (Fig. 3, Table 3).

Abnormal regional network measures in children with HCP

Among the 62 nodes included in the network, a series of nodes were found to have different regional connectivity measures with statistical significance (p < 0.05), or at least at trend level (0.05 < p < 0.1) between the controls and children with HCP in PreOp group and/or children with HCP in the PostOp group (Tables 4 and 5). Similar to the global network measurement comparisons, all the significant group differences of regional network measures between the PreOp group and the controls based on the corrected p-values using FDR approach remained statistically significant using the permutation approach, as shown in Table S2.

In the PreOp HCP group, when compared to the controls, lower nodal degree, lower betweenness centrality, higher clustering efficient, and higher local efficiency were found in the medial occipital gyrus and medial frontal gyrus bilaterally with statistical significance (p < 0.05) or at least at trend (0.05 < p < 0.1). The thalamus was also found to have lower degree and betweenness centrality bilaterally but with lower local efficiency. Higher degree and/or higher betweenness centrality were found in the cingulate gyrus, lingual gyrus, superior temporal gyrus, rectal gyrus, caudate and insular, either bilaterally or unilaterally in children in the PreOp HCP group when compared to the control group (Table 4). Among these nodes, lower clustering coefficient and/or lower local efficiency were found in the cingulate gyrus, superior temporal gyrus, and caudate in the PreOp HCP group with statistical significance or at least at trend level.

Table 4

Group comparison of regional network measures (nodal degree, betweenness centrality, clustering coefficient, and local efficiency). Between preop HCP patients and controls (all network values are residuals based on linear regression to account for age factor; all p values are FDR corrected to account for multiple comparisons across the four network measures and 62 nodes in the network; only those nodes that showed significant group difference in one or more measures are included).

Region	Degree			Betweenness centrality			Clustering coefficient			Local efficiency
	CTLn = 29	Preopn = 9		CTLn = 29	Preopn = 9		CTLn = 29	Preopn = 9		CTLn = 29	Preopn = 9
	Mean ± std	Mean ± std	t (corrected p)	Mean ± std	Mean ± std	t (corrected p)	Mean ± std	Mean ± std	t (corrected p)	Mean ± std	Mean ± std	t (corrected p)
MOG_L	0.00 ± 2.73	−5.98 ± 3.59	−5.32 (0.0002)	0.00 ± 0.03	−0.04 ± 0.01	−4.81 (0.0007)	0.00 ± 0.07	0.17 ± 0.08	6.10 (0.0001)	0.00 ± 0.05	0.11 ± 0.04	5.97 (0.0001)
MOG_R	0.00 ± 2.72	−4.32 ± 2.05	−4.38 (0.0016)	0.00 ± 0.02	−0.03 ± 0.01	−4.20 (0.0023)	0.00 ± 0.08	0.16 ± 0.12	4.61 (0.0010)	0.00 ± 0.05	0.09 ± 0.06	4.37 (0.0016)
MFG_L	0.00 ± 2.85	−5.49 ± 3.07	−4.97 (0.0005)	0.00 ± 0.01	−0.01 ± 0.01	−3.04 (0.0280)	0.00 ± 0.12	0.14 ± 0.15	3.08 (0.0273)	0.00 ± 0.06	0.07 ± 0.08	2.99 (0.0302)
MFG_R	0.00 ± 3.65	−4.05 ± 2.42	−3.11 (0.0267)	0.00 ± 0.01	−0.01 ± 0.00	−2.41 (0.0899)	0.00 ± 0.14	0.13 ± 0.17	2.31 (0.1060)	0.00 ± 0.08	0.07 ± 0.09	2.39 (0.0902)
THAL_L	0.00 ± 3.03	−8.76 ± 4.25	−6.86 (0.0001)	0.00 ± 0.02	−0.03 ± 0.02	−3.76 (0.0064)	0.00 ± 0.07	−0.04 ± 0.21	−0.88(ns)	0.00 ± 0.06	−0.16 ± 0.27	−3.09 (0.0273)
THAL_R	0.00 ± 3.60	−6.06 ± 2.65	−4.66 (0.0010)	0.00 ± 0.03	−0.03 ± 0.01	−2.81 (0.0437)	0.00 ± 0.05	0.03 ± 0.12	0.97(ns)	0.00 ± 0.06	−0.10 ± 0.19	−2.48 (0.0801)
CingG_L	0.00 ± 2.42	7.71 ± 6.44	5.44 (0.0002)	0.00 ± 0.02	0.06 ± 0.06	5.30 (0.0002)	0.00 ± 0.04	−0.09 ± 0.05	−5.08 (0.0004)	0.00 ± 0.03	−0.05 ± 0.04	−4.49 (0.0014)
CingG_R	0.00 ± 2.25	7.30 ± 4.57	6.53 (0.0001)	0.00 ± 0.02	0.05 ± 0.04	4.46 (0.0014)	0.00 ± 0.06	−0.09 ± 0.03	−4.10 (0.0028)	0.00 ± 0.04	−0.04 ± 0.02	−3.06 (0.0270)
LG_L	0.00 ± 2.02	3.90 ± 3.70	4.01 (0.0030)	0.00 ± 0.02	0.03 ± 0.03	3.53 (0.0109)	0.00 ± 0.05	−0.02 ± 0.07	−1.00(ns)	0.00 ± 0.04	−0.01 ± 0.04	−0.69(ns)
LG_R	0.00 ± 1.86	2.75 ± 3.09	3.29 (0.0193)	0.00 ± 0.02	0.01 ± 0.02	1.49 (0.3143)	0.00 ± 0.04	−0.00 ± 0.04	−0.19(ns)	0.00 ± 0.03	0.00 ± 0.03	0.33(ns)
STG_L	0.00 ± 2.55	2.50 ± 3.26	2.41 (0.0891)	0.00 ± 0.01	0.02 ± 0.02	3.01 (0.0295)	0.00 ± 0.08	−0.08 ± 0.07	−2.71 (0.0534)	0.00 ± 0.05	−0.04 ± 0.05	−2.49 (0.0810)
STG_R	0.00 ± 2.50	3.29 ± 3.38	3.16 (0.0238)	0.00 ± 0.02	0.02 ± 0.03	1.66 (0.2616)	0.00 ± 0.0	−0.04 ± 0.07	−1.77(ns)	0.00 ± 0.05	−0.02 ± 0.05	−0.88(ns)
Caud_L	0.00 ± 3.00	1.79 ± 4.68	1.36(ns)	0.00 ± 0.00	0.01 ± 0.08	2.88 (0.0380)	0.00 ± 0.19	−0.30 ± 0.21	−4.06 (0.0030)	0.00 ± 0.14	−0.21 ± 0.29	−3.08 (0.0275)
Caud_R	0.00 ± 2.20	3.50 ± 4.49	3.19 (0.0234)	0.00 ± 0.00	0.01 ± 0.01	3.95 (0.0040)	0.00 ± 0.27	−0.28 ± 0.25	−2.75 (0.0494)	0.00 ± 0.25	−0.22 ± 0.30	−2.16(ns)
Cu_L	0.00 ± 1.95	3.89 ± 3.51	4.27 (0.0020)	0.00 ± 0.01	0.02 ± 0.02	3.62 (0.0088)	0.00 ± 0.09	−0.11 ± 0.08	−3.33 (0.0177)	0.00 ± 0.05	−0.06 ± 0.04	−3.22 (0.0222)
Ins_R	0.00 ± 2.43	2.63 ± 2.13	2.90 (0.0371)	0.00 ± 0.01	0.01 ± 0.01	1.19(ns)	0.00 ± 0.11	−0.08 ± 0.10	−1.91(ns)	0.00 ± 0.07	−0.03 ± 0.06	−1.12(ns)
RG_L	0.00 ± 1.23	1.77 ± 3.28	2.43 (0.0879)	0.00 ± 0.00	0.01 ± 0.01	3.17 (0.0240)	0.00 ± 0.10	0.09 ± 0.16	1.94(ns)	0.00 ± 0.07	−0.05 ± 0.11	−1.65(ns)

In the PostOp HCP group, lower degree and lower betweenness centrality were found in the bilateral medial occipital gyrus, bilateral medial frontal gyrus, and left thalamus. Higher local efficiency and higher clustering coefficient were found in the bilateral medial occipital gyrus and medial frontal gyrus. In addition, higher degree and higher between centrality were found in the left cingulate gyrus. Lower local efficiency and lower clustering coefficient were found in the right caudate.

Discussion

Summary

In the present study graph theory analysis was applied in the assessment of structural connectivity abnormality derived from DTI network topological features in children with HCP. As compared to the controls, both HCP groups under investigation, one pre-operatively and the second post-operatively, showed significantly decreased small-worldness driven mainly by the decrease of normalized clustering coefficients. Further analysis revealed significant changes in a number of nodes in various regional network features in both groups of children with HCP.

Global network analysis

Human brain organization can be considered as a series of anatomically segregated regions, with all the regions connected by WM fibers that transfer information among them (Sporn, 2009). As described by Guye et al. (2010), there are many different types of networks based on their distinct topological network features. A regular network is a system that has all nodes connected only to their nearest neighbors with many short distance connections and without long distance connections. A random network is a system that has all its nodes connected randomly with the same probability of short and long distance connections. A “small-world” network is a system at an intermediate state between the regular network and the random network. It has many short distance connection and some long distance connections, the balance of which is quantified by the small-worldness measure. The human brain, as well as many other biological systems, has been found to possess this small-worldness feature that strikes a compromised balance between its wiring cost and resilience to pathological disturbance with the system demand for speed and synchronization.

Mathematically, the small-worldness (sigma, σ) is a ratio between two other network parameters, normalized clustering coefficient (gamma, γ) and normalized characteristic path length (lambda, λ). In theory, a decrease of small-worldness can be driven by either a decrease of normalized clustering coefficient or an increase of normalized characteristic path length, or both. Our data showed that, while the brain networks in both HCP groups preserved small-world characteristics (γ ≫ 1, λ ≈ 1, and σ ≫ 1; see γ, λ, and σ values in Table 3), the normalized clustering coefficient was significantly lower in both HCP groups, while the normalized characteristic path length was normal in the PreOp HCP group but was significantly lower in the PostOp HCP group. Although both HCP groups were found to have significantly lower small-worldness, the observations in the PreOp HCP group may suggest the beginning of brain WM injury in which the local connections have been affected as the results of the enlarged ventricles and increased ICP while the overall network level of integration has not been changed due to the redundancy at regional level. In comparison, the abnormally lower level in both normalized clustering coefficient and normalized characteristic path length in the PostOp HCP group may suggest that the integration aspect of the network has started to be affected after the surgery regardless of the diversion of CSF. However, caution should be taken in the interpretation of the differences seen in the two patient groups even if they sometimes progress in the direction that is in line with the expected recovery in response to the treatment. The results presented in this study are based on two separate cross-sectional comparisons. The different patterns in the network abnormality, e.g., the larger variance seen in the PostOp group (Fig. 3) may suggest the temporal change over time in response to surgery, but it can also be attributed to inter-subject variance at baseline. Therefore, any conclusions regarding the longitudinal change of structural connectivity over time can only be tested with a pre- and post-study design using data acquired from the same patient group.

Regional network analysis

Based on the regional network analysis in the present study, 10 brain regions (bilaterally or unilaterally) in children with HCP were found to have significant difference, or at least at trend level, in one or more of the four regional network measures, when compared to the controls. These brain regions included the media occipital gyrus, medial frontal gyrus, thalamus, cingulate gurus, superior temporal gyrus, lingual gyrus, rectal gyrus, cuneus, caudate, and insular (Tables 4 and 5). These areas are involved in many essential functions including visual and non-visual spatial processing, language, attention, working memory, learning, executive functions, motor, as well as centers for relaying various sensorimotor information for associated cortical regions (Dupont et al., 1994; Shulman et al., 1998; Renier et al., 2010). Conceivably, the abnormal level in the local connectivity measures for these brain regions may contribute to the change in normal functioning of the entire network, and thus to the change in eventual long-term outcomes. Based on the literature, the primary functional deficits in hydrocephalus patients are in the visual spatial and perceptual, visual attention, and visual motor coordination domains (Lumenta et al., 1995; Hoppe-Hirsch et al., 1998; Mataro et al., 2001; Heinsbergen et al., 2002; Persson et al., 2006, 2007). Our analysis showed that the medial occipital gyrus and thalamus (bilaterally or unilaterally) exhibited lower betweenness centrality and degree of connection in the patients. The medial occipital gyrus showed higher local efficiency and clustering coefficient in both HCP patient groups, while the thalamus showed decreased local efficiency in PreOp HCP group. The left cuneus was found to have significantly higher degree, higher betweenness centrality, lower local efficiency, and lower clustering coefficient in the PreOp group but none of these differences was significant between the PostOp HCP group and the control group. In addition, the lingual gyrus was also found to have significantly higher degree (bilateral) and higher betweenness centrality (left). It is worth noting that these regions have high potential to account for functional deficits observed frequently in hydrocephalus patients. The medial occipital gyrus is the secondary visual cortex that is involved in the visual spatial and visual perspective functioning as well as in attention and learning (Waberski et al., 2008; Tanaka et al., 2009; Schurz et al., 2013; Tu et al., 2013). The thalamus plays a key role in brain network as the gateway for relaying sensory information, including visual processing, to a wide range of cortical areas via thalamic radiation. Recent studies also demonstrate that the thalamus is an integral part of the feedback circuitry in sensory processing (see Briggs and Usery, 2008 for review). The cuneus and lingual gyrus are also known to be involved in basic visual information processing and visual memory functions (Beason-Held et al., 1998; Slotnick and Schacter, 2006; Chechlacz et al., 2012; Tamura et al., 2012; Kraft et al., 2014). It is conceivable that any disturbance to these brain regions, such as the abnormal structural connectivity measures as seen in our study, may have a ripple effect leading to disruption of normal functions in extended brain areas. Although we do not have sufficient data to establish the association between the changes in neuroanatomical substrate and the neuropsychological and behavioral outcomes in the above-mentioned domains, the data from our study based on the regional network topological features certainly point to a direction that warrant further investigation.

Compatibility and consistency of network features

It should be noted that all the network properties based on graph theoretical analysis, including global and regional measures, depend strongly on a series of factors such as parcellation scheme, imaging data acquisition protocol, and fiber tracking algorithm. As shown in a systemic review article by Zaleski et al. (2010), a network analyzed with a finer parcellation scheme (larger number of nodes) has smaller global efficiency, larger normalized clustering coefficient, slightly larger normalized path length, and larger small-worldness. Compared to conventional DTI acquisition protocol, HARDI (High Angular Resolution Diffusion Imaging), a diffusion MR imaging acquisition approach that is more sensitive to detect crossing fibers, usually generates more white matter streamlines from fiber tracking in comparison to DTI, and leads to higher global efficiency, smaller normalized clustering coefficient, smaller normalized characteristic path length, and weaker small-worldness attribute. In the present study, the brain networks from both controls and HCP patients all presented typical features of small-worldness (Watts and Strogatz, 1998; Achard et al., 2006; Humphries et al., 2006) based on the global network measures (original γ = 1.55–2.27; λ = 1.00–1.16; and σ = 1.55–2.01). After accounting for the age factor, we identified significant structural connectivity differences for the entire brain network and also in many specific brain regions. These new findings were based on the diffusion MR imaging data acquired using a 15-direction DTI protocol and the brain network with 62 cortical and subcortical regions derived from the JHU-DTI-WMPM II atlas (Oishi et al., 2008; Oishi et al., 2009; Djamanakova et al., 2013). The combination of these settings will inevitably influence the compatibilities and consistency when these results are compared with other studies.

Limitations

There are several limitations to the present study. First, the sample size for both pre- and post-op patient groups was small. Due to the typical age when CSF diversion surgeries were performed, most surgical candidates were infants. The contrast and deformation in these infants' brain images pose significant challenges for registration to a common normalized brain image template. Therefore only the older surgical candidates (>11 months) were included in the pre-operative HCP patient group. Patients at post-operative stage were older; however, many of these patients had to be excluded due to the MR image artifact caused by programmable shunt valves. Although the large deformation diffeomorphic metric mapping method used in the normalization is a highly elastic algorithm and showed good results in the present study, this step remains a challenging task in our patient population, especially when there is severe ventricle enlargement and poor image contrast. In some cases the normalization step needed to be repeated several times with manual adjustment of parameter settings. Among all the brain regions, the subcortical structures are more likely to be affected in this step due to their proximity to the enlarged ventricles in the patients. The fiber tracking and regional network measures related to these subcortical brain areas are also likely affected by the distortion. Therefore, although the present study demonstrated some interesting and promising results that were statistically significant, they all need to be replicated in a larger study population for validation. The second issue is also related to sample size, i.e., the patients in the post-operative group were scanned at either 3 or 12 month post-operative follow-up. Combining the datasets (all from different subjects) prevented us from seeing any differences between the two time points. Third, the MR sequence used in DTI data acquisition was a 15 direction sequence, which lacked the angular resolution needed for optimal white matter fiber tracking. The impact of this factor on the sensitivity and specificity of the graph analysis is unknown. However, in future study design a DTI sequence with more directions, or even a HARDI sequence should be considered. Lastly, we did not have behavioral outcome or neuropsychological evaluation data to correlate with the neuroimaging results. It should be noted that the present study is a retrospective analysis that focused on testing a new imaging analysis approach in children with hydrocephalus. The data generated from the present study should be regarded as preliminary in nature to help formulate hypotheses in future new studies.

Conclusion

The current study demonstrated that the structural connectivity analysis based on graph theory and DTI tractography is sensitive to the detection of abnormalities of brain network connectivity in children with hydrocephalus at both pre- and post-operative stages.

Table 5

Group comparison of regional network measures (nodal degree, betweenness centrality, clustering coefficient, and local efficiency) between postop HCP patients and controls. All network values are residuals based on linear regression to account for age factor; all p values are FDR corrected to account for multiple comparisons across the four network measures and 62 nodes in the network; some nodes did not show any significant group difference but are included to be consistent with Table 4.

Region	Degree			Betweenness centrality			Clustering coefficient			Local efficiency
	CTLn = 29	Postop n = 17		CTLn = 29	Postopn = 17		CTLn = 29	Postopn = 17		CTL n = 29	Postopn = 17
	Mean ± std	Mean ± std	t (corrected p)	Mean ± std	Mean ± std	t (corrected p)	Mean ± std	Mean ± std	t (corrected p)	Mean ± std	Mean ± std	t (corrected p)
MOG_L	0.00 ± 2.73	−4.06 ± 4.71	−3.71 (0.0158)	0.00 ± 0.03	−0.03 ± 0.02	−4.01 (0.0114)	0.00 ± 0.07	0.13 ± 0.13	4.42 (0.0040)	0.00 ± 0.05	0.08 ± 0.07	4.61 (0.0029)
MOG_R	0.00 ± 2.72	−2.98 ± 4.42	−2.84 (0.0543)	0.00 ± 0.02	−0.02 ± 0.02	−3.06 (0.0421)	0.00 ± 0.08	0.07 ± 0.08	2.95 (0.0482)	0.00 ± 0.05	0.04 ± 0.04	2.89 (0.0514)
MFG_L	0.00 ± 2.85	−3.76 ± 4.03	−3.71 (0.0144)	0.00 ± 0.01	−0.01 ± 0.01	−3.35 (0.0292)	0.00 ± 0.12	0.10 ± 0.13	2.60 (0.0785)	0.00 ± 0.06	0.05 ± 0.07	2.41 (0.1048)
MFG_R	0.00 ± 3.65	−2.98 ± 4.33	−2.49 (0.0906)	0.00 ± 0.01	−0.02 ± 0.02	−2.95 (0.0468)	0.00 ± 0.14	0.16 ± 0.17	3.44 (0.0247)	0.00 ± 0.08	0.08 ± 0.09	3.31 (0.0306)
THAL_L	0.00 ± 3.03	−6.74 ± 5.58	−5.32 (0.0008)	0.00 ± 0.02	−0.03 ± 0.01	−5.10 (0.0009)	0.00 ± 0.07	−0.00 ± 0.20	−0.06(ns)	0.00 ± 0.06	−0.06 ± 0.26	−1.21(ns)
THAL_R	0.00 ± 3.60	−3.16 ± 7.29	−1.97 (ns)	0.00 ± 0.03	−0.02 ± 0.02	−2.91 (0.0461)	0.00 ± 0.05	0.03 ± 0.20	0.74(ns)	0.00 ± 0.06	−0.05 ± 0.24	−0.98(ns)
CingG_L	0.00 ± 2.42	4.72 ± 7.99	2.99 (0.0475)	0.00 ± 0.02	0.03 ± 0.05	2.72 (0.0672)	0.00 ± 0.04	−0.04 ± 0.07	−2.23(ns)	0.00 ± 0.03	−0.02 ± 0.04	−1.75(ns)
CingG_R	0.00 ± 2.25	1.46 ± 6.19	1.16(ns)	0.00 ± 0.02	0.00 ± 0.04	0.44(ns)	0.00 ± 0.06	−0.01 ± 0.10	−0.23(ns)	0.00 ± 0.04	0.00 ± 0.06	0.30(ns)
LG_L	0.00 ± 2.02	−0.24 ± 4.53	−0.25(ns)	0.00 ± 0.02	−0.01 ± 0.02	−1.23(ns)	0.00 ± 0.05	0.07 ± 0.15	2.30(ns)	0.00 ± 0.04	0.04 ± 0.08	2.49 (0.09)
LG_R	0.00 ± 1.86	0.92 ± 4.22	1.02(ns)	0.00 ± 0.02	0.00 ± 0.03	0.07(ns)	0.00 ± 0.04	0.01 ± 0.11	0.45( ns)	0.00 ± 0.03	0.01 ± 0.07	0.71(ns)
STG_L	0.00 ± 2.55	−0.33 ± 2.47	−0.43(ns)	0.00 ± 0.01	0.00 ± 0.02	0.06(ns)	0.00 ± 0.08	−0.01 ± 0.08	−0.45(ns)	0.00 ± 0.05	0.00 ± 0.05	0.28(ns)
STG_R	0.00 ± 2.50	0.57 ± 3.13	0.68(ns)	0.00 ± 0.02	−0.01 ± 0.02	−1.15(ns)	0.00 ± 0.0	0.02 ± 0.06	1.28(ns)	0.00 ± 0.05	0.02 ± 0.04	1.71(ns)
Caud_L	0.00 ± 3.00	1.34 ± 5.25	1.10(ns)	0.00 ± 0.00	0.00 ± 0.01	1.54(ns)	0.00 ± 0.19	−0.14 ± 0.25	−2.22(ns)	0.00 ± 0.14	−0.11 ± 0.25	−1.84(ns)
Caud_R	0.00 ± 2.20	2.25 ± 4.76	2.19(ns)	0.00 ± 0.00	0.00 ± 0.00	2.30(ns)	0.00 ± 0.27	−0.28 ± 0.30	−3.22 (0.0348)	0.00 ± 0.25	−0.24 ± 0.37	−2.61 (0.0782)
Cu_L	0.00 ± 1.95	2.25 ± 4.84	2.23(ns)	0.00 ± 0.01	0.01 ± 0.02	2.32(ns)	0.00 ± 0.09	−0.06 ± 0.13	−1.77(ns)	0.00 ± 0.05	−0.03 ± 0.06	−2.08(ns)
Ins_R	0.00 ± 2.43	1.90 ± 3.04	2.32(ns)	0.00 ± 0.01	−0.00 ± 0.01	−0.77(ns)	0.00 ± 0.11	−0.00 ± 0.09	−0.06(ns)	0.00 ± 0.07	0.01 ± 0.04	0.67(ns)
RG_L	0.00 ± 1.23	0.49 ± 1.63	1.15(ns)	0.00 ± 0.00	0.00 ± 0.00	0.46(ns)	0.00 ± 0.10	−0.03 ± 0.16	−0.68(ns)	0.00 ± 0.07	−0.01 ± 0.07	−0.66(ns)