Altered functional and structural brain network organization in autism

1Autism spectrum disorders (ASD) are increasingly prevalent neurodevelopmental disorders characterized by atypical social behavior, including deficits in receptive and expressive language, theory of mind, and mental flexibility. Findings of functional underconnectivity between brain regions in individuals with ASD relative to matched controls have been reported as they perform a variety of cognitive tasks (see, for review). Multiple studies have found that task-independent (i.e., intrinsic) functional connectivity, including interhemispheric and default mode network (DMN) connectivity is also lower in ASD (e.g., ). Further supporting an underconnectivity theory, diffusion tensor imaging (DTI) studies have found reductions in structural white matter integrity across most major tracts (see, for review).

2In addition to reports of reduced functional connectivity within major networks (i.e., functional integration), connectivity between different networks (i.e. functional segregation) is altered in ASD. Functional brain networks become simultaneously more integrated and segregated during typical development (e.g., ) and white matter integrity increases during development (e.g., ), suggesting that brain networks in ASD may reflect ‘immature’ or aberrant developmental processes.

3Despite this array of regional and systems level findings in ASD, it is unclear how these alterations might be reflected at a network level where the brain is modeled as a network of hundreds of interacting regions composing several integrated and segregated systems. Graph theory, which describes complex systems as a set of “nodes” (i.e., brain regions) and “edges” (i.e., connections between nodes), has characterized the brain as a complex network with a hierarchical modular organization consisting of several major functional communities (i.e., visual, sensorimotor, default mode, and attentional systems; see, for review). Structural and functional brain networks exhibit robust levels of local and global efficiency (i.e., small-world properties; ) that can be quantitatively characterized using graph theoretical methods. Structural and functional graph theoretical studies have begun to map how local and global network properties change during development, aging (e.g., ) and in diseases such as schizophrenia (e.g., ) and Alzheimer's (e.g., ).

4In this study we sought to compare functional and structural connectivity in children and adolescents with ASD relative to typically developing (TD) children by characterizing local and global graph theoretical metrics of structural and functional networks using a recently validated 264-region functional parcellation scheme. We first compared simpler network connections and then characterized higher-level network properties including clustering, characteristic path length, small worldness and modularity. Additionally, since structural connectivity has been shown to correlate with functional connectivity, we wanted to determine whether structure–function correlations differed between groups and how functional and structural network properties relate to each other across development in TD and ASD individuals.

5High-functioning children and adolescents with ASD, as well as TD children and adolescents, were recruited through UCLA's Center for Autism Research and Treatment (CART) and flyers posted throughout the greater Los Angeles area. Individuals with metal implants, psychiatric or neurologic disorders, structural brain abnormalities, or known genetic conditions were excluded from participation. Informed consent and assent to participate was obtained prior to assessment according to protocols approved by the UCLA Institutional Review Board (IRB). Verbal, performance, and overall intelligence were assessed for each participant using the Wechsler Abbreviated Scale of Intelligence (WASI; ) or the full Wechsler Intelligence Scale for Children (WISC; ). High-functioning children with ASD had a prior clinical diagnosis of autism based on criteria from the Diagnostic and Statistical Manual of Mental Disorders (DSM IV), which was confirmed with the Autism Diagnostic Observation Scale (ADOS-G; ) and/or Autism Diagnostic Interview (ADI-R; ).

6A total of 60 individuals with ASD (52 males and 8 females) and 45 TD individuals (38 males and 7 females) were included in either the resting state, DTI or combined resting state/DTI final matched datasets (Table 1). After excluding subjects with excessive head motion, the resting state sample included 42 ASD subjects and 37 TD subjects and the DTI sample included 51 ASD subjects and 43 TD. Both structural and functional data were available for 35 ASD and 35 TD subjects. The three sets of matched groups did not significantly differ based on age, sex, mean/maximum head motion, or full-scale, verbal and performance IQ (Table 1).

7Twenty-two individuals with ASD and one TD individual reported the use of one or more psychotropic medications. One TD subject was using a psychostimulant. Of the subjects in our ASD sample, 12 were taking psychostimulants, 5 were taking sympatholytics, 9 were taking atypical antipsychotics, 9 were taking selective serotonin reuptake inhibitors, 3 were taking selective norepinephrine reuptake inhibitors, 3 were taking an atypical antidepressant, and 2 were taking anticonvulsants. There were no significant differences (ps > 0.30) between medicated and unmedicated ASD individuals for each of the functional and structural measures described in the following sections.

8All resting-state fMRI and DTI scans were acquired on a Siemens 3 T Trio at UCLA. A scout localizing scan was collected to help prescribe the orientation of the scans. Next, a matched bandwidth T2-weighted high-resolution echo planar scan was acquired co-planar to the functional images, which ensures identical distortion characteristics for registration purposes (Siemens 3 T Trio: TR = 5000 ms, TE = 34 ms, matrix size: 128 × 128, 19.2 cm FoV, and 36 4-mm thick slices with an in-plane voxel dimension of 1.50 × 1.50 mm). In a single session, subjects were asked to relax and keep their eyes open while a fixation cross was displayed on a white background for 6 min (T2*-weighted functional images: TR = 3000 ms, TE = 28 ms, matrix size 64 × 64, 19.2 cm FoV, and 34 4-mm thick slices (no gap), interleaved acquisition, with an in-plane voxel dimension of 3.0 × 3.0 mm). The DTI sequence consisted of 32 scans with different diffusion-weighted directions (b = 1000 s/mm2), three scans with no diffusion sensitization, at b = 0, and additional six scans at b = 50 s/mm2. Other parameters were TR = 9500 ms, TE = 87 ms, GRAPPA on, FOV = 256 mm, with 75 slices, yielding an in-plane voxel dimension of 2 × 2mm with 2-mm thick axial slices, and total scan time = 8 min 1 s.

9Imaging data were analyzed using FSL version 4.1.4 (FMRIB's Software Library, www.fmrib.ox.ac.uk/fsl; ) and AFNI (Analysis of Functional NeuroImages; ). Structural and functional images were skull-stripped using AFNI (3dskullstrip and 3dautomask). Functional volumes were motion corrected to the mean functional volume with MCFLIRT (Motion Correction using FMRIB's Linear Image Registration Tool) using a normalized correlation ratio cost function and sinc interpolation. Translations and rotations in the x, y, and z dimensions were calculated from volume to volume and averaged to generate mean and max relative displacement values, which did not significantly differ between the final matched groups (Table 1). Subjects with a single displacement (combined translational and rotational movements) greater than 2.5 mm (13 ASD and 5 TD) were excluded prior to further analyses and not included in the final samples. Images were spatially smoothed using a Gaussian kernel of FWHM 5 mm. A band pass filter (0.1 Hz > t > 0.01 Hz) was applied to the data in order to minimize the effects of cardiac and respiratory fluctuations. The 6 rigid body motion parameters and average white matter (WM), cerebrospinal fluid (CSF), and global time-series and their temporal derivatives were then regressed out of the data. The WM and CSF time-series reflected signal from subject-specific regions of interest created using FAST (FSL's Automatic Segmentation Tool). Given the recent concerns regarding the effect of motion in resting state fMRI connectivity, in addition to matching the groups by mean and maximum relative head motion, we also regressed out individual volumes with large signal intensity changes (i.e., motion spikes) by creating additional nuisance regressors that modeled individual time points with greater than half of a standard deviation change in global signal intensity.

10One major methodological hurdle in graph theory approaches to neuroimaging concerns how to define the nodes of the network. Most studies have used anatomical atlases (e.g., ) or individual voxels (e.g., ) as nodes. However, anatomical atlases include relatively large regions that are likely to contain multiple functional regions, which can distort/obscure true properties of the network by mixing distinct signals. Conversely, voxel-wise parcellation approaches can be biased by artificially strong local connections. A whole-brain parcellation scheme was recently created based on a large meta-analysis of fMRI studies combined with whole brain functional connectivity mapping. This set of 264 putative functional regions was shown to more accurately represent the information present in the network (i.e., it was better at detecting previously characterized functional networks such as dorsal and ventral attention subnetworks) relative to voxelwise and atlas-based parcellation approaches. Therefore, we chose this set of 264 regions for whole-brain parcellation. For each subject, 5-mm radius spheres based on the MNI coordinates of these 264 regions were registered to functional space (12 DOF, affine, and correlation ratio cost function) through registration from the MNI 152 template to the high-resolution echo-planar (12 DOF, affine, and mutual information cost function) using FSL's Linear Image Registration Tool (FLIRT). We then correlated timeseries between each of the 264 brain regions and z-transformed correlation coefficients in order to generate 264 × 264 whole brain functional connectivity matrices for each subject. Graph theoretical metrics and statistics were computed with Matlab (The Mathworks, Natick, MA) using the Brain Connectivity Toolbox.

11Before comparing graph theoretical properties, we sought to investigate pairwise differences in connection strengths as a function of the network's modularity, which refers to the set of subnetworks or distinct communities that exist within the network as a whole. Constructing the most representative modularity partition is an active area of research (http://arxiv.org/abs/1206.4358); however, at present there is no single, agreed-upon method for choosing the most representative partition. Here we used the Louvain modularity algorithm applied to the unthresholded functional connectivity matrix averaged across all subjects after removing all negative weights. Although there is some recent work regarding how to incorporate negative weights into graph theoretical metrics, we chose to use well-established algorithms that use only positive connections in the calculation of graph theoretical metrics. We chose a random partition containing 4 modules (the most common number of modules identified over 100 runs of the algorithm) as our representative partition and reorganized the order of nodes in the functional connectivity matrix by this modular organization for visualization purposes (Fig. 1A). This representative partition was also used to determine whether each connection was a within- or between-module connection for additional calculations described in the next paragraph. The similarity of this chosen modularity partition with 99 other modularity iterations of the group average matrix was calculated using normalized mutual information (NMI; ). Additionally, given the controversy regarding the accuracy of comparing modularity of group average matrices with individual subject matrices, we compared the similarity of each individual's modular organization with that of our representative group average modular organization with NMI.

12In order to interrogate pairwise differences in connection strengths between the two groups, two-sample t-tests were performed for every z-transformed connection strength value (Fig. 1C). The significance was set at p < 0.05, uncorrected, for these initial exploratory analyses. If there was significantly lower connectivity in a connection for the ASD or TD group (vs. the other group) that had an average correlation value below zero, it was categorized as stronger negative connectivity for that group (as opposed to stronger positive connectivity for the other group, which is mathematically equivalent) given that negative or anticorrelations are likely to represent real phenomenon and stronger connectivity for one group simply represent weaker negative connectivity for the other group. The number of connections differing between groups was assessed for each of the identified modules both for within-module positive connections and between-module negative connections (Fig. 1D). Numbers of within-module positive and between-module negative connections differing between groups were compared and displayed as a function of the connection's average z-transformed correlation value (Fig. 2A) and Euclidean distance between regions (Fig. 2B).

13As there is no rationale for using a particular cutoff for functional connectivity strength to determine whether an edge exists in a functional network, we compared local and global network properties over a range of functional connection thresholds. Thresholding a network based on correlation strength can yield different network sparsities (number of existing edges divided by number of possible edges), which influence network properties and can bias a comparison of graph metrics between groups. In fact, we found that at higher z-correlation thresholds the TD group had a higher average sparsity (Fig. 3A). Therefore, we chose to equalize network sparsity between subjects by taking an equivalent percentage of the strongest positive connections (negative connections were ignored) for each subject and binarizing the network weights before calculating graph theoretical metrics. Binarization is a common step in functional graphs (e.g., ) in order to preserve only the strongest (most probable) functional connections and treat these connections equivalently. We examined functional network properties between 15% and 32% sparsity. The upper threshold of 32% was chosen because the weakest edge it includes corresponds to a correlation coefficient of 0.15, which is the minimum correlation needed to be statistically significant (p < 0.05) across 120 functional images. At lower sparsity levels, network properties begin to break down as the network becomes fragmented. Therefore, we chose 15% sparsity (corresponding to a minimum correlation coefficient of 0.34 (p <.001)) as the low end of the range based on the requirement that all individual subject graphs be fully connected (Fig. 3B).

14We focused on 6 global graph theoretical metrics (see for formulas of these metrics). These metrics were: clustering coefficient (CC), which measures how much neighbors of a node are connected to each other and is closely related to local efficiency; characteristic path length (CPL), which is the average number of edges needed to get from any node in the network to any other node in the network and is inversely related to global efficiency; normalized CC and CPL (lambda and gamma), which are calculated as the ratios of CC or CPL to the average CC or CPL from simulated random networks; small worldness, which is the ratio of lambda to gamma; and modularity Q values, which represent the proportion of within-module edges in the network minus within-module edges calculated from a similar random network. For the calculation of lambda and gamma, we randomized networks by starting with a true network and then performing random double edge swaps with the constraint that these swaps must maintain the connectedness of the network. This algorithm preserved the degree of each node in the true network and was performed with the randmio_und_connected.m script in the Brain Connectivity Toolbox. One hundred of these random networks were calculated for each subject. Lambda and gamma were calculated using the mean of the C and L from the random networks. Since modularity Q values can vary based on random differences in module assignments from run to run, Q values were averaged over 100 iterations of the algorithm. All metrics were averaged across 15% to 32% sparsities in 1% increments to generate average values for each metric given the smooth curve across the sparsity range (Fig. 3). Two sample t-tests were performed on these metrics between subjects at each sparsity level (Fig. 3C–H) and for metrics averaged across sparsity levels (Table 2). To correct for multiple comparisons across the 6 metrics, False Discovery Rate (FDR q < 0.05; ) was applied. For each node, clustering coefficients, participation coefficients and betweenness centrality were also averaged across sparsity levels for each subject and compared between groups (Fig. 4A). Betweenness centrality measures how often the shortest path goes through a given node while participation coefficients reflect how much a node interacts with nodes in different communities and each roughly corresponds to global metrics of characteristic path length and modularity, respectively. Differences in nodal metrics are shown at more stringent (FDR: q < 0.05) and less stringent thresholds (p < 0.05, uncorrected).

15Individual volumes with gross motion artifacts were excluded from further analysis and subjects with excessive motion (greater than 8 volumes (20%) with motion artifacts) were not included in final samples (6 ASD and 3 TD). Motion and eddy current correction was performed on the diffusion-weighted images using eddy_correct in FMRIB's Diffusion Toolbox (FDT), while MCFLIRT was used to quantify mean and maximum relative motion (Table 1), which did not differ between groups. Dtifit was used to fit a diffusion tensor model to the data at each voxel and calculate voxelwise Fractional Anisotropy (FA) values for each subject. Whole brain deterministic tractography was then performed using the fiber assignment by continuous tracking (FACT) algorithm in Diffusion Toolkit (http://trackvis.org/dtk). We sought to boost the likelihood of detecting longer fibers between spatially separate spherical ROIs by relaxing constraints on our tractography algorithm. Therefore, tractography was carried out by propagating fibers from each voxel with a maximum turn angle of 50° and without an FA cutoff. The spatial separation of the ROIs effectively acts as a filter and offsets the reduced constraints placed on tractography as, with greater distance, it becomes less likely that spurious fibers will continue to propagate and connect distant ROIs. Fibers were smoothed using a spline filter. Fibers shorter than 5 mm were excluded as this corresponds to 2 voxels, for which a turn angle cannot be determined.

16We used the same set of 264 coordinates from to generate 10 mm radius spheres in MNI space. Dilating spheres to 10 mm radii (relative to 5 mm radii spheres for functional nodes) ensured inclusion of nearby white matter fibers given that nodal coordinates were centered in gray matter. This set of nodes covers 50.6% of all white matter voxels based on FSL's white matter tissue priors thresholded at 50%. Additionally, on average 60.9% of the voxels in each ROI were white matter voxels. These 264 masks were transformed to each subject's diffusion space (12 DOF, affine, and correlation ratio cost function) through registration to the hires image (12 DOF, affine, and mutual information cost function). In order to generate edges between nodes of structural networks, the number of fibers connecting each region was counted. A fiber was defined as connecting two regions if one fiber endpoint terminated within one region and the other endpoint terminated within the other region. This process was repeated using all 264 regions as seeds in order to derive a 264 × 264 whole brain structural connectivity matrix for each subject, using custom software written for this purpose (UCLA Multimodal Connectivity Package; http://github.com/jbrown81/umcp). Additionally, average FA and mean diffusivity (MD) were calculated for each connection.

17The Louvain modularity algorithm was run on the group average unthresholded fiber connectivity matrix. The order of nodes in the fiber connectivity matrix was reorganized based on a representative modularity partition with 9 modules (Fig. 5A). The similarity of this representative average modularity partition was compared with each individual matrix's modularity partition as well as 99 additional runs of the modularity algorithm on the average matrix using normalized mutual information.

18Two sample t-tests (p < 0.05, uncorrected for initial exploratory analyses) were performed on fiber counts, FA, and MD for every connection after masking by connections that have an average of 5 or more fibers (5.75% of all possible connections; Fig. 5B,C) in order to minimize false positive connections. Connections differing between groups for number of fibers and MD were compared as a function of average fiber count and Euclidean distance.

19For structural networks, we examined the same six global network properties as functional networks (CC, CPL, lambda, gamma, small worldness and modularity Q values) averaged between 5% and 8.5% sparsity in 0.5% increments. Structural networks were then binarized in order to maintain maximum comparability to equivalent functional networks. A sparsity level of 5% represented the minimum sparsity level at which every subject's graph was fully connected (Fig. 6A) and 8.5% represented the average unthresholded sparsity of all subject's structural matrices (Fig. 6B). Two sample t-tests were performed on these six metrics between subjects for averaged metrics (with FDR correction) as well as at each sparsity level.

20Fiber counts of every connection with an average of at least 5 fibers were correlated with functional connectivity strengths for each of the 35 ASD and 35 TD subjects (Fig. 7). Additionally, fiber count/functional connectivity correlations were computed for within- and between-module connections and specifically for within-module connections with lower levels of functional connectivity as identified in Fig. 1C. These structure–function correlations were z transformed, then compared between groups (with two-sample t-tests).

21We ran an exploratory principal component analysis (PCA) on the six average functional global graph metrics and the six average structural graph global metrics across all 70 subjects with PASW Statistics 18, Release Version 18.0.3 (SPSS, Inc., Chicago, IL). Values for the first four of 12 total components were computed for each subject and compared between groups (two sample t-tests, with FDR correction), and correlated with chronological age after regressing out mean relative motion (with FDR correction; Table 3; Fig. 8A,B). The two components that significantly differed between groups were also tested for correlation with symptom severity (as measured by the social and communication subscales of the ADOS and ADI; with FDR correction) within the ASD group after regressing out mean relative motion and age (Table 3; Fig. 8C–D).

22Renderings were generated from scripts in the UCLA Multimodal Connectivity Package (http://github.com/jbrown81/umcp) and through the UCLA Multimodal Connectivity Database (http://umcd.humanconnectomeproject.org), which use matplotlib (http://matplotlib.sourceforge.net) and networkX (http://networkx.lanl.gov).

23All of the connectivity matrices used in this study are freely available for download at the UCLA Multimodal Connectivity Database (; http://umcd.humanconnectomeproject.org).

24Over the course of 100 runs of the modularity algorithm on the average functional connectivity matrix, four communities were detected over 90% of the time. The normalized mutual information between the representative partition and 99 additional iterations of the modularity algorithm was 0.73 +/− 0.05 (mean +/− standard deviation; NMI ranges from 0 to 1, with 1 representing a perfect similarity). There was moderate similarity between the representative group average modularity partition and modularity partitions of individual subjects (ASD group NMI = 0.285 +/− 0.07 and TD group NMI = 0.288 +/− 0.08; p = 0.86). The average NMI between randomized individual matrices and the representative group average modularity partition was equal 0.04, suggesting that the similarity between individual matrices and the group average was above chance.

25The order of nodes in the matrix was reorganized to reflect the community structure of the representative group average modularity partition (Fig. 1A). The four communities corresponded to visual, sensorimotor and default systems as well as a largely frontal system corresponding to the task positive control/attention network (color boxes in Fig. 1A and displayed in 3D brain space in Fig. 1B).

26We first examined pairwise differences in the connectivity matrices by directly comparing correlation strengths between groups for each connection (Fig. 1C) and separating differences based on within- and between-community connections. We found that the TD group exhibited 5.4 times as many stronger (p < 0.05, uncorrected) within-module positive connections as the ASD group (Fig. 1C,D). This was most pronounced in the default (265 (10.1%) connections stronger for TD group vs. 15 (0.5%) stronger for ASD), visual (107 (7.5%) connections stronger for TD group vs. 7 (0.5%) stronger for ASD) and sensorimotor systems (84 (2.5%) connections stronger for TD group vs. 33 (1.0%) stronger for ASD; Fig. 1D, left). There were a similar number of stronger within-module connections for the attention/control network (34 (2.3%) connections stronger for TD group vs. 41 (2.3%) stronger for ASD). Additionally, the TD group exhibited 4.4 times as many stronger (p < 0.05, uncorrected) negative (i.e., weaker) between-module connections. This was most prominent for connections between other systems and the default (670 (4.8%) for TD > ASD vs. 152 (1.1%) for ASD > TD) system, but was also true for visual (383 (3.4%) for TD > ASD vs. 135 (1.2%) for ASD > TD), sensorimotor (479 (3.2%) for TD > ASD vs. 154 (1.0%) for ASD > TD), and attention (276 (2.4%) for TD > ASD vs. 165 (1.5%) for ASD > TD) systems (Fig. 1D, right). Thus, there was a pattern of weaker within-network positive connectivity and weaker between-network negative connectivity for children and adolescents with ASD.

27We sorted within and between-module differences as a function of average correlation strengths (Fig. 2A). Connections where the TD group had stronger positive within-module connectivity tended to have higher average correlation strengths than connections where the ASD group had stronger within-module connectivity (TD = 0.26 +/− 0.19, ASD = 0.16 +/− 0.17, p = 0.0002). Between-module connections where the TD group had stronger negative connectivity were more negative than the connections where the ASD group had stronger negative connections (TD = − 0.16 +/− 0.09, ASD = − 0.08 +/− 0.06, p < 0.0001). We found no significant differences (all p > 0.25) for the average Euclidean distance of connections that differed between groups for stronger positive within-module connectivity or stronger negative between-module connectivity (Fig. 2B).

28There were group differences in nearly all graph theoretical metrics for functional networks over a range of network sparsities (Fig. 3) and averaged across sparsity levels (Table 2). Clustering coefficient was significantly lower (FDR q < 0.05) in the ASD group (Fig. 3C) and although lambda was lower in the ASD group at higher sparsity levels (Fig. 3D), there was only a trend for lower average gamma. Both CPL and gamma were lower in the ASD group over the entire range of sparsities (FDR q < 0.05; Fig. 3E,F) and averaged across sparsity levels (Table 2). Both TD and ASD subjects had functional networks in the small world range (the ratio of lambda to gamma being greater than 1.2). However, small worldness was not significantly different between groups. Modularity (Q values) was significantly lower (FDR q < 0.05) in the ASD group at every sparsity level and averaged across sparsity levels (Fig. 3H). In addition to averaging metrics across this sparsity range, metrics were also integrated across sparsity levels as in, and then compared between groups. There was a 0.99 correlation between averaged and integrated metrics. Additionally, there were no alterations in any of the results reported above when using integrated metrics.

29Given the significant between-group differences in global metrics for CC, CPL and modularity, we sought to determine whether specific nodes contributed to these global differences. Therefore, we compared nodal measures of local interconnectivity, hubness, and connection diversity by calculating each node's clustering coefficient, betweenness centrality, and participation coefficient between groups averaged over the same range of thresholds (Fig. 4). We report the number of nodes with significant between-group differences (p < 0.05 uncorrected and FDR corrected; q < 0.05, p < 0.0013). ASD subjects had lower nodal CC in 21 visual (4 FDR: right occipital fusiform gyrus and left and right inferior lateral occipital cortex), 20 default (3 FDR: medial prefrontal cortex, ventromedial prefrontal cortex and left frontal orbital cortex), and 10 sensorimotor nodes (1 FDR: left superior parietal lobule; Fig. 4A). Participation coefficients were higher for the ASD group in 26 default (3 FDR: medial prefrontal cortex and left frontal orbital cortex), 10 sensorimotor (3 FDR: left postcentral gyrus, left superior parietal lobule and brainstem) and 9 attention (0 FDR) nodes (Fig. 4B). There were no differences in nodal betweenness centrality that survived FDR correction.

30The Louvain modularity algorithm detected between 8 and 10 communities for the average fiber connectivity matrix over 100 runs. Nine communities were detected in over 80% of the runs and these communities corresponded to sets of lateralized nearby brain regions (Fig. 5B). The average fiber structural connectivity matrix for all TD and ASD subjects is shown in Fig. 5A, after reordering the nodes by the community structure of a representative modularity partition.

31The calculated normalized mutual information between the representative structural modularity partition and 99 additional iterations of the modularity algorithm was 0.84 +/− 0.04. There was also high similarity between the representative group average modularity partition and modularity partitions of individual subjects (ASD group NMI = 0.66 +/− 0.05; TD group NMI = 0.68 +/− 0.05; p = 0.14).

32We first examined the connectivity matrices by directly comparing the number of fibers, average FA, and average MD values for each connection between groups after masking for regions that contained an average of at least 5 fibers (corresponding to 5.75% of all possible connections). We found that the ASD group had 4.2 times as many connections with significantly (p < 0.05, uncorrected) more fibers than the TD group (106 ASD > TD vs. 25 TD > ASD; Fig. 5C). We also found that the ASD group had 1.6 times as many connections with lower FA (67 TD > ASD vs. 41 ASD > TD) and 6.2 times as many connections with higher MD (112 ASD > TD vs. 18 ASD > TD; Fig. 5D).

33The average number of fibers or Euclidean distance of the connection did not differ for connections where the ASD group had more fibers compared to connections where the TD group had more fibers (number of fibers: TD > ASD = 30.0 +/− 26.6, ASD > TD = 26.3 +/− 25.9, p = 0.53; Euclidean distance: TD > ASD = 25.9 +/− 12.8, ASD > TD = 34.8 +/− 26.3, p = 0.11). Connections where the TD group had higher white matter integrity (lower MD) had a higher average number of fibers than connections where the ASD group had higher white matter integrity (ASD > TD = 32.5 +/− 28.8, TD > ASD = 11.8 +/− 6.0, p = 0.003), but did not differ based on Euclidean distance (ASD > TD = 30.3 +/− 21.5, TD > ASD = 35.2 +/− 7.5, p = 0.34).

34Although gamma (normalized characteristic path length) was similar for structural and functional networks (~ 1.2 for structural vs. ~ 1.1 for functional), lambda (normalized clustering coefficient) was much higher in structural networks (~ 5.4 for structural vs. ~ 2.2 for functional). Therefore, structural networks displayed higher levels of small worldness compared to functional networks in both TD and ASD groups. Measures for average structural CC, gamma, CPL, lambda and small worldness did not significantly differ between groups (Table 2). Modularity Q values were higher in the TD group on average, but this did not survive FDR correction.

35Given the previous reports of modularity decreasing with age and global efficiency increasing with age in structural networks, we ran post-hoc analyses correlating these metrics with chronological age in each group. Higher modularity in the TD group was actually driven by the younger TD participants, whereby, controlling for motion, modularity was significantly negatively correlated with age in the TD group (r = − 0.41, p = 0.008; Fig. 6D) yet was only trending toward a negative correlation with age in the ASD group (r = − 0.24, p = 0.08) although the interaction was not significant (p = 0.37). Similarly, there were no group differences for CPL or lambda, but age was negatively correlated with CPL and lambda in the TD group (controlling for motion: CPL: r = − 0.34, p = 0.03; lambda r = − 0.31, p = 0.04), and CPL and lambda were positively correlated with age in the ASD group (controlling for motion: CPL: r = 0.22, p = 0.12; gamma: r = 0.16, p = 0.25; Fig. 6C) whereby there was a significant group by age interaction for CPL (p = 0.01) and gamma (p = 0.02).

36When comparing correlations between fiber counts and functional connectivity strength between groups, we found that both groups exhibited moderate, yet highly significant (all subjects p < 0.001), levels of structural–functional connectivity correlations (TD: r = 0.32 +/− 0.03, ASD: r = 0.32 +/− 0.04, p = 0.77 for the group difference). Furthermore, there were no group differences when structure–function correlations were assessed for both within- and between-functional module connections (within-module: TD: r = 0.28 +/− 0.04, ASD: r = 0.28 +/− 0.05, p = 0.98 for the group difference and between-module TD: r = 0.28 +/− 0.05, ASD: r = 0.26 +/− 0.05, p = 0.32 for the group difference) or specifically for within-module connections exhibiting lower levels of functional connectivity (TD: r = 0.30 +/− 0.12, ASD: r = 0.31 +/− 0.11, p = 0.77 for the group difference).

37To identify key factors underlying correlated graph metrics and to better understand relationships between structural and functional network properties, we entered the six functional and six structural average global graph metrics for all 70 subjects into an exploratory principal component analysis. We only examined the first four components, as they explained the vast majority (88%) of the variance in the data (Table 3). The first component (accounting for 33.9% of the variance) broadly weighted functional metrics positively and structural metrics negatively. This component was significantly lower in the ASD group (covarying for mean head motion, p = 0.009, FDR q < 0.05) and negatively related to symptom severity, as measured by the ADI social subscale (covarying for age and mean head motion, r = − 0.4, p = 0.01; FDR q < 0.05; Fig. 8C). The first component was also positively correlated with age in both groups (covarying for mean motion; all: r = 0.24, p = 0.04, TD: r = 0.28, p = 0.11 and ASD: r = 0.30, p = 0.08; Fig. 8A). The second component weighted all functional and structural metrics positively, and although it did not differ between groups, there was a significant interaction with age (covarying for mean head motion; p = 0.02), whereby the second component was significantly positively correlated with age in the ASD group (r = 0.35, p = 0.04) and slightly negatively correlated with age in the TD group (r = − 0.24, p = 0.16). The third component, positively weighting functional CC/modularity and negatively weighting functional CPL, did not differ between groups. The fourth component, positively weighting structural modularity and negatively weighting structural CPL, was significantly lower in the ASD group (p = 0.007; FDR q < 0.05) and was negatively correlated with symptom severity as measured by the ADOS social and communication subscales (covarying for age and mean head motion, ADOS social: r = − 0.46, p = 0.005; FDR q < 0.05, Fig. 8D).

38Previous neuroimaging studies on ASD have reported reduced functional and structural connectivity both within and between specialized brain systems, suggesting ASD is a network disorder. Here we expand upon previous findings of lower functional and structural connectivity in ASD by characterizing higher-level network properties using tools derived from the physics of complex networks. We report alterations in community organization of functional networks, as well as in the balance of local and global efficiency within and between structural and functional networks in children and adolescents with ASD relative to their typically-developing counterparts.

39We detected robust reductions in positive functional connectivity within major functional systems (i.e., functional integration) in individuals with ASD. Reduced functional connectivity was most prominent in the default system, consistent with multiple studies that have found reduced DMN connectivity in ASD. However, we also found weaker connectivity within visual (largely secondary areas) and sensorimotor systems, supporting more widespread alterations in functional connectivity as found by,, and. Relatively few alterations were observed in the frontal attention/cognitive control network, which might reflect relatively intact cognitive skills in high-functioning individuals with ASD.

40Interestingly, individuals with ASD also show reduced negative (i.e., more positive) connectivity between systems. Consistent with previous findings in the task positive and default mode networks in ASD, weaker negative connectivity between communities suggests that specific functional systems are less distinct or functionally segregated from one another. Although there is some controversy regarding the proper interpretation of negatively correlated brain regions when using global signal regression (GSR; ), anticorrelations are detected without GSR and GSR maximizes the specificity of positive resting-state correlations in real and simulated data. Interestingly, reduced negative connectivity was recently shown to be useful for diagnostic classification of autism in analyses without GSR. Therefore, although it is unclear whether widespread differences in negatively connected regions are exaggerated by GSR, differences in negative connectivity between distinct functional systems are likely important for understanding ASD neurobiology.

41Although most previous functional connectivity studies of ASD have reported underconnectivity of long-range (i.e., anterior–posterior or interhemispheric) connections, it has also been widely hypothesized that ASD may be related to overconnectivity of short-range connections. Previous neuroimaging studies have found increased short-range connections in neurotypical children versus adults but findings are somewhat mixed in individuals with ASD. Consistent with a recent study we found that even short-range functional connections are reduced in ASD. Of course, this does not exclude the possibility that local connections at the neuronal or minicolumnar level could be enhanced in ASD.

42In examining graph metrics of functional networks, we found that individuals with ASD had lower clustering (i.e., local efficiency), especially in nodes within the default systems and secondary visual areas. Individuals with ASD displayed a less robust modular organization (i.e., communities were less distinct) and there was a tendency for nodes in the default and sensorimotor systems to interact more with other communities as measured by higher nodal participation coefficients. Finally, we found that functional brain networks in individuals with autism had shorter average path lengths (i.e., higher levels of global efficiency) as well as normalized characteristic path lengths. Randomly connected networks tend to have short path lengths suggesting the possibility that higher global efficiency in functional networks may simply reflect a less organized or more random distribution of functional edges. This is consistent with a study finding decreased complexity or increased randomness in resting-state fMRI timeseries of individuals with ASD.

43Previous functional graph theory studies in typical development did not find differences in local or global efficiency between children and adults. However, in addition to decreases in long-range connectivity, these developmental studies reported increased local connectivity, which may explain these null findings. Additionally, the extent to which these previously reported developmental differences are attributable to motion artifacts is unclear given that subtle motion spikes tend to reduce long range connectivity yet increase local connectivity. Given our careful consideration of head motion through regression of motion spikes and covarying for motion at the group level, as well as the fact that we found both reduced long- and short-range connectivity in ASD, it is unlikely that our findings are related to between-group differences in motion. Although future work is needed to further examine developmental changes in the context of more stringent motion correction, our findings in ASD are somewhat consistent with the studies in typical development reporting reduced integration and segregation of functional systems in children relative to adults. Therefore, although functional networks in ASD may be ‘immature’ in some ways (i.e., reflect an earlier developmental stage as far as reduced integration/segregation of major systems), they may also be fundamentally different from neurotypical individuals from a network perspective (i.e., reduced local efficiency yet increased global efficiency).

44For structural connectivity measures derived from diffusion MRI, we found reduced integrity in short- and long-range white matter tracts in ASD in line with previous studies (e.g., ). We found more robust differences in MD than FA, which has been reported in several previous DTI studies. Interestingly, despite the fact that white matter integrity was generally reduced, we found evidence for increased fiber counts in ASD, which may relate to early reports of increased regional white matter and more recent reports of increased fiber counts in certain tracts in ASD. Although white matter integrity is lower in children compared to adults, fiber counts increase during development. Therefore, like functional networks, some alterations in ASD may reflect immaturity, while other alterations are likely to reflect aberrant processes.

45Structural networks displayed high levels of local and global efficiency in both the TD and ASD groups. Given previous reports of decreasing modularity and increasing global efficiency of structural networks with development, we examined the relationship between age and modularity/global efficiency in each group. We found that in the TD group, modularity sharply decreased with age whereas global efficiency increased with age, consistent with previous reports. In the ASD group, modularity decreased at a slower rate and, contrary to findings in the TD group, global efficiency actually decreased with age. It should be noted that global efficiency in structural networks likely reflects a different underlying substrate than global efficiency in functional networks given the physical wiring costs of structural networks. Thus, despite similar levels of local and global efficiency in structural networks across both groups, it appears as though network efficiency does not appropriately shift from a more local to a more distributed pattern during development in individuals with ASD.

46One potential limitation of our analyses is that we chose to use a set of spherical, functionally-based ROIs instead of more traditional block or atlas based regions, which would have allowed for the inclusion of deeper white matter and, relatedly, increased signal to noise in the tractography analyses. This decision was made in order to allow for a more direct comparison between structural and functional connectivity. Additionally, previous structural connectivity studies have used fibers that terminate at this gray–white boundary because they are the most reliable/likely estimates of cortico-cortical connectivity. However, future work should incorporate other parcellation schemes that include deeper white matter while also allowing for a direct comparison of structure and function.

47When relating structural and functional connectivity, we found that measures of fiber counts and functional connectivity strength were moderately positively correlated in both groups with no group differences regardless of whether the connections were within or between modules or whether we only included connections with lower levels of functional connectivity. This finding, in addition to the fact that we generally saw higher fiber counts in ASD, suggests that alterations in functional connectivity in ASD are not directly related to alterations in fiber organization.

48In order to relate structural and functional network properties, we performed a principal component analysis. Interestingly, we found that the largest underlying factor inversely weighted structural and functional network properties. This component inversely weighted local and global efficiency (i.e., positively weighted both CC and CPL) within functional and structural networks and was positively correlated with age. Although preliminary, this finding suggests that structural networks become more globally efficient, yet less locally efficient, during development while functional networks display a relative inverse pattern. This component was reduced in ASD and inversely related to social and communicative behavior, suggesting that the balance between structural and functional network properties is related to social impairments in ASD. Further highlighting differential age-related trajectories for functional and structural network properties, the second component, which positively weighted both structural and functional metrics, decreased with age in the TD group while it increased with age in the ASD group. Interestingly, a previous study of multiple sclerosis found that disruption of white matter pathways actually leads to increased functional connectivity in multiple networks including the DMN, which further highlights a divergence between the structural and functional connectomes. However, it should be noted that here we found inherently different sparsity ranges for structural and functional networks. Direct comparison of structural to functional graphs at different sparsities is problematic given that graph theoretical metrics can vary as a function of network sparsity. Our PCA method attempts to overcome the confound of direct comparison between differentially sparse structural and functional graphs in order to integrate information across modalities.

49Finally, the fourth component, which positively weighted local and global efficiency in structural networks, was reduced in ASD and inversely related to social and communicative symptom severity. Therefore, an underlying factor positively influencing both local and global efficiency in structural networks may also relate to disrupted social behavior in ASD.

50Future studies should characterize younger and/or lower functioning individuals with ASD since our findings are limited to high-functioning children and adolescents with ASD. For example, studies examining infants at high risk for ASD may be useful for developing biomarkers to aid in earlier diagnosis and treatment. Future studies may also benefit from advances in imaging acquisition, more flexible modeling approaches, and large-scale studies involving collaboration between institutions. Additionally, comparisons with other neuropsychiatric disorders, and teasing apart underlying mechanisms such as genetic risk factors will all be crucial for a more complete characterization of brain network abnormalities in ASD.

51To our knowledge, this is the first study to use complex network analyses to examine both structural and functional brain networks in autism. We found significant reductions in local efficiency and modularity within several functional networks. ASD children and adolescents also displayed atypical age-related changes in the balance of local and global efficiency between structural and functional networks. Further, this imbalance was related to the severity of socio-communicative deficits in individuals with ASD. Our findings suggest that complex network modeling of structural and functional brain organization will yield a better understanding of the neural basis of ASD and other neuropsychiatric disorders. Ultimately, a more cohesive framework for understanding brain alterations in ASD may inform the design of more sophisticated diagnostic tools and targeted interventions.