Disrupted anatomic networks in the 22q11.2 deletion syndrome.

The 22q11.2 deletion syndrome (22q11DS) is an uncommon genetic disorder with an increased risk of psychosis. Although the neural substrates of psychosis and schizophrenia are not well understood, aberrations in cortical networks represent intriguing potential mechanisms. Investigations of anatomic networks within 22q11DS are sparse. We investigated group differences in anatomic network structure in 48 individuals with 22q11DS and 370 typically developing controls by analyzing covariance patterns in cortical thickness among 68 regions of interest using graph theoretical models. Subjects with 22q11DS had less robust geographic organization relative to the control group, particularly in the occipital and parietal lobes. Multiple global graph theoretical statistics were decreased in 22q11DS. These results are consistent with prior studies demonstrating decreased connectivity in 22q11DS using other neuroimaging methodologies.

Introduction

Recent advances in molecular genetics provide an intriguing set of putative biomarkers associated with increased liability to numerous psychiatric diseases (Kendler, 2013). For schizophrenia specifically, the genetic variant of greatest effect is a well-established copy number variant at 22q11.2 (International Schizophrenia Consortium, 2008). A hemizygous deletion at this locus, typically of about 1.5-3 Mb (Shaikh et al., 2000), results in a broad spectrum of craniofacial, cardiac, endocrine, neurologic, and psychiatric manifestations now commonly known as the 22q11 deletion syndrome or 22q11DS (Shprintzen, 2008). People with 22q11DS have increased risk for several psychiatric conditions including attention-deficient hyperactivity disorders, mood disorders, anxiety, and autism-spectrum disorders (Tang et al., 2013). Perhaps most strikingly, persons with 22q11DS have an estimated 50% lifetime prevalence of subthreshold psychotic features and a 25-fold increased risk of psychotic-spectrum disorders relative to the general population (Bassett and Chow, 1999, Murphy et al., 1999).

Idiopathic schizophrenia is increasingly considered a disorder of altered brain connectivity (Alexander-Bloch et al., 2010). Prior studies using diffusion tensor imaging (DTI) suggest lower fractional anisotropy and higher diffusivity in individuals with schizophrenia relative to typically developing controls (Kubicki et al., 2007). DTI studies in 22q11DS are less common, but also generally suggest lower anisotropy relative to typically developing controls (Barnea-Goraly et al., 2003, Jalbrzikowski et al., 2014, Simon et al., 2005, Villalon-Reina et al., 2013). Resting state fMRI studies also suggest altered connectivity in 22q11DS (Debbané et al., 2012).

The analysis of covariance patterns within structural data using graph theoretical models provides an additional method to explore the multivariate substrates of complex neuropsychiatric diseases (Alexander-Bloch et al., 2013, Bassett et al., 2008). By observing differential correlation patterns among anatomic regions of interest, the underlying neural network can be inferred (Chen et al., 2012, He et al., 2007). Network modeling with cortical thickness data, in particular, reproduces known patterns of axonal connectivity (Lerch et al., 2006). In the current study, we compare network structure in 22q11DS to a large group of typically developing controls via a multivariate analysis of cortical thickness. We hypothesized that several graph theoretical statistics would be decreased in 22q11DS relative to typically developing controls, mirroring prior brain network analyses using other imaging modalities.

Methods

Sample

The 22q11DS sample was drawn from a prospective study, Brain-Behavior and Genetic Studies of the 22q11DS at University of Pennsylvania and Children's Hospital of Philadelphia (CHOP). Participants were recruited from the “22q and You Center” at CHOP and through social media. Inclusion criteria were: age ≥ 8, English proficiency, estimated IQ > 70 by clinical testing and the Wide Range Achievement Test IV (Wilkinson and Robertson, 2006), and stable medical status. Exclusion criteria were: pervasive developmental disorder or IQ < 70, medical disorders that may affect brain function (e.g., uncontrolled seizures, head trauma, CNS tumor and infection) or visual disability (e.g., blindness); appropriate subjects older than 12 years of age were considered for neuroimaging. The 22q11.2DS neuroimaging subsample consisted of 48 individuals (mean age 20.1 years ± SD 4.5, 52% male) and is similar to a sample described previously (Schmitt et al., 2014a). More than half of the 22q11DS (54%) had significant psychosis spectrum symptoms, including 48% with sub-threshold (positive, negative or disorganized) symptoms and the remaining with threshold psychotic disorders, including 2 with schizophrenia/schizoaffective disorder. Psychiatric disorders in general were common; 27% of subjects had a history of major depression, 38% with an anxiety disorder, and 35% with ADHD. Only seven (15%) subjects had no significant psychosis spectrum symptoms or lifetime history of other psychiatric disorders.

Our typically developing control group was obtained from the Philadelphia Neurodevelopmental Cohort (PNC), a prospective sample of children and young adults aged 8–21 years, recruited through CHOP. Details on the neuroimaging sample are described elsewhere (Satterthwaite et al., 2013). A total of 370 typically developing subjects without psychiatric diagnoses or any significant psychopathology were included (mean age 14.7, sd 4.1, 48% male).

Image Acquisition

High-resolution axial T1 weighted magnetization prepared rapid acquisition gradient (MPRAGE) echo was acquired, with the following parameters; TR/TE 1810/3.51 ms; TI 1100 ms; FOV 180 × 240 mm; effective resolution 1 mm3. 22q11DS and PNC participants were all scanned using the same 3 Tesla MRI scanner (TIM Trio; Siemens, Erlangen, Germany), pulse sequence parameters and 32-channel head coil. Board certified technologists in the Department of Radiology at the Hospital of the University of Pennsylvania performed and read all scans.

Image processing

In order to obtain measures of cortical thickness, raw data were imported into FreeSurfer version 5.0 (http://surfer.mgh.harvard.edu). FreeSurfer's surface-based image processing pipeline is described extensively elsewhere (Dale et al., 1999, Fischl and Dale, 2000, Fischl, 2012, Fischl et al., 1999). Briefly, for each subject, image intensity was normalized to account for magnetic field inhomogeneity. The skull and other non-brain tissues were removed (Ségonne et al., 2004). Preliminary segmentation was performed using a connected components algorithm. The surface boundary was then covered with a polygonal tessellation and smoothed, resulting in high-resolution vertices over both cerebral hemispheres. A deformable surface algorithm was employed to identify the pial surface. The cortical surface model was reviewed and manually edited if necessary. Cortical parcellation was performed based on a combination of local curvature information, the probability of a regional label at a given location in surface-based atlas space, and contextual information (Fischl et al., 2004) based on the Desikan atlas (Desikan et al., 2006). This parcellation scheme resulted in 68 total regions of interest (ROIs). For each ROI, cortical thickness was calculated by averaging the distance between the pial surface and the gray/white boundary.

Statistical Analysis

Data were then imported into the statistical programming environment R (R Core Development Team, 2012). For each ROI, linear and nonlinear effects of age, sex, race, age-sex interactions, and average global cortical thickness were controlled via regression. Correlation matrices were constructed from the residuals for each group separately and compared using the Jennrich and Mantel tests (Jennrich, 1970, Mantel, 1967). Correlational patterns were visualized using the heatmap.2 function (gplots package) with superimposed hierarchical cluster analysis using Euclidean distances (Warnes et al., 2015). Hierarchical clustering requires no a priori specification of the number of clusters present in the data, but rather orders relationships based on a distance function (Hastie et al., 2011). Heatmap.2 performs agglomerative clustering, which represents a stepwise bottom-up strategy that recursively groups the most related structures until a single cluster remains. In addition to reorganizing the data such that ROIs with similar correlational patterns are spatially proximal in the matrix, a dendrogram also is produced that shows the level of similarity among the ROIs; the shorter the path along the dendrogram between two ROIs, the more similar their patterns of correlations.

Undirected graph theoretical models were then constructed for 22q11DS and TD groups separately. The use of binary networks decreases variability in network statistics and is particularly advantageous for relatively small sample sizes (Cheng et al., 2012). Since density of a network can influence many network statistics (Bassett et al., 2008), we examined the network properties of both groups separately using a range of correlational thresholds to define significant edges, similar to methods described in Liu et al. (2008); the analysis pipeline is summarized in Fig. 1. Undirected graphs were constructed for each group separately in R using the igraph package (Csárdi and Nepusz, 2006, Csárdi and Nepusz, 2015). Several global network statistics (described below) were calculated for each threshold separately using existing functions in the igraph or qgraph packages (Epskamp et al., 2012). Standard deviations for each statistic-threshold combination were calculated for each group separately via bootstrap with 1000 replicates. In order to provide a summary statistic independent of individual threshold selection, area under the curve (AUC) analysis was then performed for each network statistic (Long et al., 2013). AUC was mathematically defined as:where n is the number of thresholds, Y(Tk) represents the value of a network statistic at threshold k, and I is the interval between thresholds. Standard deviations in AUC were calculated via bootstrap with 1000 replicates, and significant group differences in AUC identified via t-test. Additionally, since two networks at a fixed correlation threshold may not have identical densities, we repeated these analysis over a range of fixed sparcities.

Fig. 1

Summary of the statistical pipeline. Correlation matrices for corrected measures of cortical thickness (A) were constructed for TD and 22q11DS groups separately. Significant edges were then identified either via serial thresholding or using the PCIT algorithm (B). Finally, undirected networks for each group were generated in igraph (C). The same pipeline was used for subsequent bootstrap and permutation analyses.

As an alternate strategy, significant edges were identified using the PCIT algorithm (Reverter and Chan, 2008), which combines information theory and first order partial correlation coefficients in order to minimize spurious edges. The algorithm was originally designed to identify meaningful gene co-expression networks, but is mathematically generalizable to most correlational data and is available in the R package ‘pcit’ (Watson-Haigh et al., 2010). In order to ensure that differences in density did not bias network statistics (Bassett et al., 2008), the density of the TD group was set to match the 22q11DS group by adjusting the correlation threshold such that graph density was equal between groups (Bassett et al., 2008, He et al., 2009). In practice, following application of the PCIT algorithm the number of edges was nearly equivalent (213 in TD group and 214 in 22q11DS), even prior to adjustment. Undirected graphs were then constructed for each group separately using the igraph package.

Using these graph models, we then explored global network architecture through several measures of connectivity (mean betweenness, modularity, clustering coefficient, average path length, and smallworldness) commonly used in graph theory (Rubinov and Sporns, 2010, van Wijk et al., 2010). For a given graph G with V vertices and E edges, these statistics are defined as follows:

Mean Betweenness

The shortest path length between two vertices (i.e. ROIs) vi and vj represents the number of connecting edges for the minimum path between them. The betweenness centrality of a vertex v is defined as the proportion of shortest paths between other vertex-vertex pairs that traverse it (Brandes, 2001, Freeman, 1978). Thus, it quantifies how a node influences the connectivity between other nodes. For a vertex v along a path connecting vertices i and j, the betweenness can be expressed mathematically as:where σ is the number of shortest paths between i and j. The mean betweenness is simply the average betweennesses for all vertices.

Modularity

The modularity Q is defined as the proportion of edges that are located within network submodules in excess of that expected by random chance (Clauset et al., 2004, Girvan et al., 2002). Thus, it is a measure of the community structure of a network, with higher values representing more densely connected vertices within modules. For a graph with binary adjacency matrix A and m total edges and for all pairwise combinations of vertices i and j:where k represents the degree (i.e. number of connections) of vertices i and j. δ = 1 if i and j are in the same module and 0 otherwise.

Clustering Coefficient

The global clustering coefficient (C, also known as transitivity) represents the average proportion of a vertex's neighbors that are also connected to one another (Barrat et al., 2004, Watts and Strogatz, 1998). The clustering coefficient of a single vertex represents the ratio of closed triplets divided by the total number of all triplets (open and closed); in other words, it represents the proportion of all possible connections with neighboring vertices that are actually realized. Strogatz’ classic example is that within a social network, C represents the fraction of one's friends that are also friends with each other. The clustering coefficient can be generalized to the network level. For all vertices with degree > 2:

Average Path Length

The average (or characteristic) path length is simply the average shortest path length between all pairwise combinations of nodes (Watts and Strogatz, 1998). For shortest path length l(vi,vj) and n vertices, the average path length lμ is defined as:

Small Worldness

Small world networks are characterized by dense local interconnectivity and short local path lengths between two vertices (Humphries and Gurney, 2008, Watts and Strogatz, 1998). Multiple prior studies have suggested that the human brain has intrinsic small world properties (Bullmore et al., 2009) that may be at least in part genetically-mediated (Schmitt et al., 2008). The smallworldness S of a graph is defined as the ratio of its clustering coefficient C to its average minimum path length L (i.e. the average number of edges between vertices), normalized to a randomly generated graph with an identical number of vertices and edges:

Smallworldness was calculated using the ‘qgraph’ package in R (Epskamp et al., 2012). Group differences in global network statistics were estimated with difference scores. In order to test for group differences from PCIT analyses, we performed permutation analyses to estimate the sampling distribution of network statistics. This was accomplished by randomly assigning individuals to two groups (equal in proportion to the original data) and repeating the analysis pipeline described above, calculating a difference score between the two groups (Fisher, 1935, Nichols and Holmes, 2003). The null distribution was estimated using 10,000 replicates.

Subgroup Analyses

Given the age difference between the 22q11DS group and our younger controls, we repeated our analysis using subsets of both groups that were better age matched by only including control subjects older than 13 and 22q11DS subjects < 24. A total of 37 individuals with 22q11DS and 150 TD controls were included. There were no significant differences in mean age between this control subgroup (mean age 18.2 years ± SD 2.26) and the group with 22q11DS (18.7 ± 3.19).

Results

When comparing correlation matrices as a whole, there were statistically significant differences in correlational patterns between groups (Jennrich χ2 = 4422.4, p-value < 0.0001; Mantel p-value < 0.0001). Four distinct clusters were identified in the TD group: 1) a cluster including most ROIs in the frontal lobe, 2) a parietal-occipital cluster (with subclusters for occipital and parietal lobe ROIs), 3) a temporal lobe-insular cluster, and 4) a cluster including the ROIs of the cingulate and precentral gyrus (Fig. 2). Rather than observing hemispheric segregation, ROIs were usually tightly correlated with their contralateral analogs. Overall cross-trait correlations between ROIs in 22q11DS were weaker; clustering patterns also were less striking, with less distinct associations between ROIs in spatial proximity.

Fig. 2

Results of hierarchical cluster analysis for cortical thickness for TD and 22q11DS groups. Positive correlations are shown in green, negative correlations in red. In the TD group, ROIs were strongly clustered by lobar anatomy: 1) frontal, 2) parieto-occipital, 3) insulo-temporal, and 4) limbic. Cross-trait correlations in 22q11DS were weaker, as was the degree of anatomic clustering.

Graph analysis showed reduced mean betweenness, modularity, clustering coefficient, average path length, and smallworldness in 22q11DS relative to the TD group over a range of correlation thresholds (Fig. 3) and sparcities (Fig. 4). AUC analysis confirmed overall reduced network measures in 22q11DS relative to the control group (Fig. 3, Fig. 4, insets), with statistically significant group differences in AUC (p < 0.0001).

Fig. 3

Group differences in global cortical thickness network properties as a function of correlational threshold. Network statistics are shown for typically developing (black circles) and 22q11DS (red triangles) separately with error bars representing standard deviations. Insets display results of AUC analysis for each network statistic; there were significant group differences (p < 0.0001) for all statistics.

Fig. 4

Group differences in global cortical thickness network properties as a function of network sparcicity. Network statistics are shown for typically developing (black circles) and 22q11DS (red triangles) separately with error bars representing standard deviations. Insets display results of AUC analysis for each network statistic; there were significant group differences (p < 0.0001) for all statistics.

Network statistics for PCIT-determined graphs are shown in Table 1. When these graphs were visualized, the 22q11DS group was subjectively less tightly connected (Fig. 5A). In TD subjects, anatomic regions generally segregated by anatomic lobe, with particularly tight interconnections within and between the occipital lobes. Subjectively, the TD group had higher modularity relative to 22q11DS. In contrast, local connectivity patterns were less striking in 22q11DS, with reduced clustering in the occipital lobe, parietal lobe, and limbic structures in particular. Quantitative metrics of connectivity confirmed these observations, with significantly lower global measures of modularity, betweenness, clustering, and smallworldness (p < 0.001, Fig. 3B).

Table 1

Comparison of Network Statistics for graphs with significant edges defined by the PCIT algorithm.

	Full Sample		Age-Matched
TD	22q11DS	TD	22q11DS
Mean Betweeness	79.54	61.15	76.47	57.44
Modularity	0.5792	0.4967	0.5756	0.4269
Smallworldness	3.484	2.606	3.284	2.739
Clustering Coefficient	0.4558	0.3147	0.4128	0.2680
Average Path Length	3.374	2.825	3.273	2.714

Fig. 5

Disrupted connectivity and modularity in 22q11DS. Graph models of cortical thickness networks for TD and 22q11DS (top). Nodes represent 68 regions of interest color-coded by lobar anatomy. Node shape indicates laterality (square = right, circle = left). Global connectivity statistics (bottom). Black dots represent difference scores (TD - 22q11DS) for measures of network cohesion. The null distribution was estimated empirically via permutation, with mean (open dot) and 95% confidence intervals.

Subgroup analysis

Repeat analysis using age-matched subgroups produced largely similar findings to the full dataset (Supplementary Results). Again, there were statistically significant group differences in correlational patterns (Jennrich test p < 0.0001; Mantel test p < 0.0001). Correlation maps subjectively appeared more tightly clustered in the TD group relative to 22q11DS. AUC for network was significantly increased (p < 0.0001) in the TD group for smallworldness, betweenness, and modularity relative to 22q11DS; clustering coefficient was higher in 22q11DS in the subgroup model, although this difference was not statistically significant (p = 0.0629). The PCIT algorithm produced similar findings in the subgroup analyses with significantly increased network statistics for all network measures in the TD group relative to 22q11DS.

Discussion

The present study provides further evidence that the 22q11.2 deletion results in global disruptions in cerebral anatomic connectivity relative to the typically developing population. To our knowledge, this represents the first exploration of network structure in 22q11DS using covariance patterns in cortical thickness data. In addition to identifying statistically significant reductions in global structure, we observed particular decreases in the strength of parieto-occipital and limbic networks compared to typically developing controls.

Our data largely support prior observations using other modalities (Debbané et al., 2012, Ottet et al., 2013b, Padula et al., 2015, Scariati et al., 2016, Simon et al., 2005). For example, Ottet et al. reported significant reductions in occipito-occipital and limbic connectivity in 22q11DS using DTI tractography in 30 individuals with 22q11DS and 30 typically developing controls, as well as a 10% reduction in the total number of fibers (Ottet et al., 2013a). Ottet et al. later expanded these analysis using graph theoretical models in a sample of 46 participants with 22q11DS and 48 matched controls, identifying significantly increased global path lengths and reduced global efficiency (Ottet et al., 2013b). Using resting state fMRI, Debbané et al. reported that individuals with 22q11DS had weaker connectivity in visuospatial, frontotemporal, and sensorimotor networks, as well as within the default mode network (Debbané et al., 2012). Alterations to the dorsal stream have long been postulated to explain the observed deficits in visuospatial and mathematical ability in 22q11DS, and several prior structural studies have identified reduced parieto-occipital volumes in people carrying the deletion (Eliez et al., 2000, Jalbrzikowski et al., 2013, Kates et al., 2001, Schmitt et al., 2014a). Studies on cortical thickness specifically have shown particular reductions in 22q11DS within parieto-occipital structures, the region whose correlational patterns appear most disrupted in our analyses (Bearden et al., 2007, Jalbrzikowski et al., 2013, Schmitt et al., 2014a). For example, Bearden et al. compared cortical thickness in 21 individuals with 22q11DS to 13 typically developing controls and found the largest region of cortical thinning in the occipital pole extending into the superior parietal lobe.

Recent graph theoretical studies on idiopathic schizophrenia using resting state fMRI, DTI, and anatomic networks also have generally suggested aberrant network connectivity (Alexander-Bloch et al., 2013, Van Den Heuvel and Fornito, 2014). For example, prior resting state fMRI studies have reported reductions in smallworldness, average clustering, and global efficiency (Lynall et al., 2010), reductions in clustering and smallworldness (Liu et al., 2008), and reductions in local efficiency, clustering, and smallworldness (Alexander-Bloch et al., 2010) in people with schizophrenia. Graph theoretical studies using DTI have shown loss of distributed connectivity within the parieto-occipital network and its connections to the frontal lobe through the cingulum and cingulate cortex (van den Heuvel et al., 2010, Zalesky et al., 2011). Employing methods similar to ours, Zhang et al. used Freesurfer-derived measures of cortical thickness to explore cerebral connectivity differences between 101 subjects with schizophrenia and 101 matched controls (Zhang et al., 2012). In addition to identifying global alterations in smallworldness, they reported relative reductions in betweenness in schizophrenia, localizing to the bilateral midline parietal lobes, parahippocampal/lingual gyri, right superior frontal gyrus, and operculum. Notably, multiple studies have suggested that people with 22q11DS have structural anomalies including most of these regions, with particular involvement of parasagittal cortical structures (Jalbrzikowski et al., 2013, Schaer et al., 2008, Schmitt et al., 2014a, Schmitt et al., 2014b).

Conclusions

Graph theoretical analysis of cortical thickness in 22q11DS supports prior findings of disrupted structural and functional connectivity in this condition, with many similarities to network architecture in schizophrenia. These findings provide further evidence that specific anatomic anomalies in 22q11DS should not be considered in isolation, but rather as genetically mediated alterations to global neurodevelopmental patterning.

Limitations

The current study has several limitations that should be considered. First, in order to maximize the available sample, there were group differences in age that were controlled statistically. Although we attempted to minimize the contributions of age by controlling for both nonlinear effects and interactions, some bias cannot be entirely excluded. Our subgroup analysis with better age matching supports the general conclusions obtained from the full sample. Second, there are substantial differences in IQ between typically developing individuals and individuals with 22q11DS that are difficult to control for. A study design including a second control group matched for IQ could improve specificity although in practice this design is difficult to implement. Additionally, IQ-matching may have its own limitations since it could result in selecting under-achieving controls and over-achieving patients (Resnick, 1992). Third, the use of group statistical correlational patterns limits our ability to assess for individual differences. Fourth, multivariate anatomic approaches infer connectivity patterns rather than visualizing them directly. However, the use of anatomic data does have some relative advantages compared to other modalities, notably the availability of larger sample sizes and increased precision of measurement, and therefore should be considered complementary to other multivariate approaches.