Building functional connectivity neuromarkers of behavioral self-regulation across children with and without Autism Spectrum Disorder.

Behavioral self-regulation develops rapidly during childhood and struggles in this area can have lifelong negative outcomes. Challenges with self-regulation are common to several neurodevelopmental conditions, including Autism Spectrum Disorder (ASD). Little is known about the neural expression of behavioral regulation in children with and without neurodevelopmental conditions. We examined whole-brain brain functional correlations (FC) and behavioral regulation through connectome predictive modelling (CPM). CPM is a data-driven protocol for developing predictive models of brain-behavior relationships and assessing their potential as 'neuromarkers' using cross-validation. The data stems from the ABIDE II and comprises 276 children with and without ASD (8-13 years). We identified networks whose FC predicted individual differences in behavioral regulation. These network models predicted novel individuals' inhibition and shifting from FC data in both a leave-one-out, and split halves, cross-validation. We observed commonalities and differences, with inhibition relying on more posterior networks, shifting relying on more anterior networks, and both involving regions of the DMN. Our findings substantially add to our knowledge on the neural expressions of inhibition and shifting across children with and without a neurodevelopmental condition. Given the numerous behavioral issues that can be quantified dimensionally, refinement of whole-brain neuromarker techniques may prove useful in the future.

Introduction

Childhood is a period when critical skills such as behavioral self-regulation are acquired and continually refined (Dajani and Uddin, 2015). Behavioral regulation skills allow children to succeed at school as well as socially (Diamond, 2013; Faja and Nelson Darling, 2018; Blair and Raver, 2015). At school, these skills enable children to follow teacher instructions despite distractions and to modulate arousal in a novel and demanding environment, which both contribute directly to their academic performance (Graziano et al., 2007; Blair, 2002). Socially, a child who has better self-regulation will have less frequent and less intense shifts in behavior, leading to more positive relationships with teachers and peers (Blair and Raver, 2015; Graziano and Hart, 2016). Relatively weaker behavioral self-regulation on the other hand associates with greater daily-life challenges and an increased risk for psychiatric diagnoses (Faja and Nelson Darling, 2018; Leung et al., 2016; de Vries and Geurts, 2015).

Children with neurodevelopmental conditions (NDCs) are known to struggle with behavioral self-regulation and exhibit problems with executive functions more generally (Schmitz et al., 2006; Schmitt et al., 2018; Ting and Weiss, 2017; Shaw et al., 2014). Executive function difficulties have been both theoretically and experimentally linked to the diagnostic symptoms of NDCs, such as most strikingly in Autism Spectrum Disorder (ASD): difficulty with normed social interactions and communication; circumscribed interests and repetitive behaviors; and both internal (anxiety, depression) and external (aggression, hyperactivity) emotional struggles (Faja and Nelson Darling, 2018; Leung et al., 2016; Ting and Weiss, 2017; Turner, 1999; Mosconi et al., 2009; de Vries and Geurts, 2012; Lopez et al., 2005; Yerys et al., 2009; D’Cruz et al., 2013; Miller et al., 2015; Lieb and Bohnert, 2017). Yet, little is known about the brain correlates of behavioral regulation across typically developing children and children with NDCs.

Behavioral regulation skills begin to develop in early childhood, rapidly increase during elementary school years and continue to improve through adolescence and adulthood (Anderson, 2002; Dick, 2014; Cepeda et al., 2001). Children may have trouble in multiple aspects of behavioral regulation subdomains, including inhibitory control processes, cognitive flexibility, and emotional control processes (Diamond, 2013; Faja and Nelson Darling, 2018; Chan et al., 2008; Gioia et al., 2000). Inhibitory control is the ability to suppress interfering distractions and prepotent motor responses (Diamond, 2013; Faja and Nelson Darling, 2018; Nigg, 2000). Cognitive flexibility, which is typically measured using set-shifting, refers to the readiness with which one can switch from one task or mindset to another (Diamond, 2013; Miyake et al., 2000; Armbruster et al., 2012). Finally, emotional control, also termed ‘cognitive control of emotion’ or emotion regulation, is the processes by which we influence which emotions we have, when we have them, and how we experience and express them (Ochsner et al., 2012; Gross, 2002). Studies indicate that reduced inhibition (Faja and Nelson Darling, 2018; Schmitt et al., 2018; Geurts et al., 2014; Voorhies et al., 2018), reduced cognitive flexibility (Van Eylen et al., 2015; Geurts et al., 2009) and reduced emotion control (Ting and Weiss, 2017; Richey et al., 2015; Berkovits et al., 2017) may all be issues for children with ASD.

The brain networks implicated in these three subdomains of behavioral regulation are heavily intertwined. Common pathways - largely identified through task-based fMRI studies - most prominently revolve around areas in the prefrontal cortex (PFC), such as the ventrolateral PFC (vlPFC), the dorsolateral PFC (dlPFC), ventromedial and dorsomedial PFC (vmPFC/dmPFC), as well as the anterior cingulate cortex (ACC). For example, data suggests that the vlPFC supports reflexive reorienting, motor inhibition, and action updating (Levy and Wagner, 2011), the selection of the most efficacious response set when confronted with a task requiring various possible responses (Badre and Wagner, 2006; Dippel and Beste, 2015), and deliberately increasing or decreasing negative affect (Ochsner et al., 2004). Enhanced dlPFC activity and also enhanced activity in medial prefrontal structures is often observed in conflict paradigms, which require inhibitory control of prepotent but incorrect responses and set-shifting to reframe the problem (Oldrati et al., 2016); such paradigms include emotional conflict paradigms where the preprotent responses that need to be controlled are emotional (Rohr et al., 2016; Etkin et al., 2010; Egner et al., 2008). In addition, areas in the parietal and temporal lobes, such as the inferior parietal lobule (IPL), the superior temporal gyrus (STS) and the temporal pole, as well as limbic structures, most notably the insula and amygdala, are known to play important roles in cognitive flexibility (Niendam et al., 2012) and emotion control (Rohr et al., 2016; Zaki et al., 2012; Rohr et al., 2015; Ferri et al., 2016).

Despite a dominant focus on prefrontal brain regions in the behavioral regulation literature, distributed regions and their interactions have also been implicated. For example, changes in inter-regional synchrony, i.e. functional connectivity or correlation (FC), during emotional control processes have been observed between IPL and vmPFC, and STS and dmPFC (Rohr et al., 2016), as well as between parietal structures and the amygdala (Ferri et al., 2016). Meta-analytic evidence (Niendam et al., 2012) lends further support to the involvement of a distributed network including parietal and temporal areas and the importance of FC. Importantly, FC can be measured both during behavioral regulation tasks, as in the literature noted above, as well as at rest, where it is thought to reflect intrinsic functional brain organization that is reflective of an individual’s attributes, including cognitive and behavioral traits (Rohr et al., 2015; Raichle, 2015; Rohr et al., 2013). FC work on this topic in NDCs in sparse; in ASD, distributed network patterns have been associated with social symptoms (Lake et al., 2019), but less is known about distributed FC in relation to behavioral regulation.

It should be noted that while the neural substrates of executive functioning and behavioral regulation have been extensively studied in neurotypical adults (Rohr et al., 2016, 2015) and adolescents/young adults with affective or developmental disorders (Etkin et al., 2010; Kana et al., 2007; Solomon et al., 2009), less is known generally about the neural expression of behavioral regulation in children and children with NDCs such as ASD (Dajani and Uddin, 2015; Yerys et al., 2009): previous work was conducted in relatively small samples (total N range = 24–38) and with varying success (Ambrosino et al., 2014; Shafritz et al., 2015; Lee et al., 2009). The neural signatures of many NDCs, and in particular ASD, have been elusive due to heterogeneity in these conditions, so taking a dimensional approach to look at specific features may be more promising (Lake et al., 2019; Ameis et al., 2016; Uddin et al., 2017). In particular, whole-brain FC signatures that associate with aspects of behavioral regulation are not well studied. This is despite the widespread repercussions of suboptimal behavioral regulation that persist into adulthood, and the enormous potential of whole-brain FC profiles as ‘neuromarkers’ for diagnosis and individually tailored treatment.

Here we use connectome predictive modelling (CPM), a data-driven protocol for developing predictive models of brain–behavior relationships from FC data using cross-validation (Shen et al., 2017), to examine whole-brain linear associations of behavioral regulation and their utility as ‘neuromarkers’. Using data from children with and without ASD, who together show a range of behavioral regulation skills, we hypothesize (a) that FC models of behavioral regulation can be built across an aggregate sample including data from children with and without an NDC, (b) that frontoparietal, limbic and default mode networks underlie behavioral regulation, and (c) that neuromarkers built in a subset of the sample can then be used to predict behavioral regulation scores in another, unseen subset of children.

Methods

Participants

For this study we used two datasets from the Autism Brain Imaging Data Exchange II (ABIDE-II) database (Di Martino et al., 2017). These sites were chosen because they included Behavior Rating Inventory of Executive Function (BRIEF (Gioia et al., 2000)) behavioral regulation scores and resting state fMRI data from children aged 8–13 years with and without ASD. Specifically, we analyzed the datasets collected at Georgetown University (GU) and the Kennedy Krieger Institute (KK), which are publicly available under fcon_1000.projects.nitrc.org/indi/abide/abide_II.html. At both sites, children were introduced to an MRI simulator first and given the opportunity to familiarize themselves with the experience of undergoing an MRI scan with their eyes open. Potential participants were excluded if they had a history of neurological or psychiatric disorders (the latter only in the typically developing (TD) group), contraindications for MRI, or if they had other medical problems that prevented participation. Nine participants did not have behavioral data on file and were excluded from the analysis. Participants’ data were further evaluated for outliers in behavioral measures, high inconsistency scores on the BRIEF or excessive motion on the fMRI scans. For the behavioral measures, outliers were defined as > 3 SD from the mean. One participant was excluded due to this criterion in the subdomain of emotion control and shifting, and 6 in the subdomain of inhibition. An additional participant was excluded due to this criterion in analyses that involved the social responsiveness scale (SRS). Five participants were excluded because of inconsistency scores >7. Twenty-six participants had excessive head motion on their fMRI scan (>4 mm maximum absolute displacement). The final samples and participant characteristics are given in Table 1.

Table 1

Participant characteristics in the GU site, the KK site and the combined sample. Means and standard deviations (in brackets) are given for the total samples comprised of both TD and ASD participants, which were used to build the models, as well as for TD and ASD participants separately. Motion (mm) refers to the absolute maximum displacement at any timepoint in the resting-state fMRI scan prior to motion mitigation and denoising procedures. n = number of participants; m = male; f = female; L = left-handed; A = ambidextrous; R = right-handed; IQ = Intelligence Quotient; FIQ = full scale IQ; VIQ = verbal IQ; PIQ = performance IQ; SRS = Social Responsiveness Scale (total score). SRS, inhibition, shifting and emotion control are given as T scores. * denotes deviating numbers in the Inhibition models, for which an additional six outliers (>3 SD in score; all ASD participants) were removed. ** SRS scores were not available for some participants in the KK site (140 out of 145 TD and 47 out of 49 ASD); 1 participant with ASD from the GU site was removed as an outlier in SRS scores. † denotes a significant difference between TD and ASD children (see Supplementary Table S1 for p-values).

	GU Site			KK Site			Combined Sample
	Total	TD	ASD	Total	TD	ASD	Total	TD	ASD
n	82	46	36	194	145	49	276	191	85
Age	10.8 (1.6)†	10.5 (1.7)	11.2 (1.4)	10.3 (1.3)	10.3 (1.2)	10.3 (1.5)	10.5 (1.4)	10.4 (1.3)	10.7 (1.5)
Sex (m/f)	54/28†	23/23	31/5	128/66	93/52	35/14	182/94†	116/75	66/19
Handedness (L/A/R)	9/0/73	3/0/43	6/0/30	12/16/166	10/10/125	2/6/41	21/16/239	12/10/168	8/6/71
Motion	1.1 (0.7)	1.1 (0.7)	1 (0.7)	1.4 (0.9)	1.2 (0.8)	1.8 (0.9)	1.3 (0.9)†	1.2 (0.8)	1.5 (0.9)
FIQ	120.2 (14.2)	121.7 (14.2)	118.3 (14.3)	111.9 (12.7)†	114.6 (10.4)	103.7 (15.5)	114.3 (13.7)†	116.3 (11.8)	109.8 (16.6)
VIQ	121.8 (14.9)	121.8 (15.5)	121.7 (14.4)	116 (13.7)†	118.2 (11.8)	109.4 (16.9)	117.7 (14.3)†	119.1 (12.8)	114.6 (16.8)
PIQ	115.9 (13.7)	116.9 (13.4)	114.1 (14.4)	109.4 (12.8)†	110.8 (12.1)	105.4 (14.3)	111.1 (13.4)†	112.2 (12.6)	108.4 (14.8)
SRS**	56.3 (17.2)†	43.9 (6.4)	72.6 (12.5)	51.5 (16.1)†	43.2 (5.3)	76.4 (10.7)	53 (16.5)†	43.4 (5.6)	74.8 (11.5)
Inhibition*	50.1 (9.8)†	44.9 (6.3)	56.8 (9.5)	47.9 (10)†	43.8 (6.1)	61.7 (8.4)	48.6 (10)†	44.1 (6.1)	59.4 (9.2)
Shifting	53.8 (14.6)†	44.6 (6.8)	65.5 (13.5)	49.7 (13.9)†	43.4 (7.1)	68.5 (12)	50.9 (14.2)†	43.7 (7)	67.2 (12.7)
Emotion Control	51.5 (11.7)†	45.2 (7.7)	59.6 (11.1)	47.5 (11)†	43. 5 (7.1)	59.4 (11.8)	48.7 (11.3)†	43.9 (7.3)	59.5 (11.5)

Cognitive and behavioral assessment

The cognitive and behavioral data used for analysis included IQ, handedness, ASD symptoms and behavioral regulation scores. IQ was measured using the Wechsler Intelligence Scale for Children – Fourth or Fifth Edition (WISC IV (Wechsler, 2003); WISC-V (Wechsler, 2014)) across the KK site, while in the GU site it was assessed using the WISC-IV or the Wechsler Abbreviated Scale of Intelligence (WASI (Wechsler, 1999)). Handedness was determined with the Edinburgh Handedness Inventory (Oldfield, 1971) in the KK site; in the GU site handedness was obtained by self- and parent-report. As an estimate of ASD symptoms, social responsiveness was measured using the Social Responsiveness Scale (SRS), Edition 1, version 1 (Constantino and Gruber, 2007) in the GU site and either in Edition 1, version 1 or Edition 2, version 1 (Constantino and Gruber, 2012) in the KK site. Behavioral regulation in the datasets was assessed with the BRIEF (Gioia et al., 2000), a parent assessment of executive function behaviors for children aged 5–18 years. BRIEF subscales provide measures of three domains of behavioral regulation, which are labelled “inhibit”, “shift”, and “emotional control”. The “inhibit” subscale assesses the ability to resist impulses and to stop one’s own behavior” (sample item: “acts wilder or sillier than others in groups (birthday parties, recess)”). The “shift” subscale assesses the ability to move freely from one situation, activity, or problem to another; to tolerate change, and to switch or alternate attention (sample item: “resists or has trouble accepting a different way to solve a problem with schoolwork, friends, chores, etc.”). Finally, the “emotional control” subscale assesses the ability to regulate emotional responses appropriately (sample item: “overreacts to small problems”). T-scores for both the SRS and the BRIEF were used for analysis.

Analysis of cognitive and behavioral measures

To assess differences in characteristics (demographics, cognitive abilities and outcomes of interest) between children with and without ASD, as well as between groups within and across the two sites, t-tests were computed. Further, to assess the relationship between behavioral regulation scores, ASD symptom scores and IQ, Pearson's correlations were computed. One-way ANOVAs were used to assess potential differences in behavioral regulation scores, ASD symptom scores and IQ in relation to handedness. Behavioral analyses were carried out using SPSS 22 (Chicago, IL).

MRI data acquisition parameters

Data were acquired on a 3 T Siemens Magnetom TrioTim at the GU site, and on a 3 T Philips Achieva at the KK site. Children were instructed to keep their eyes open at both sites. Functional images were acquired using a gradient-echo EPI sequence in 43 axial slices (154 volumes, TR = 2000 ms, TE = 30 ms, FA = 90, matrix size 64 × 64, voxel size 3 × 3 × 3 mm³; duration: 5.14 min) at the GU site, and in 47 axial slices (128 volumes, TR = 2500 ms, TE = 30 ms, FA = 75, matrix size 96 × 96, voxel size 2.67 × 2.67 × 3 mm³; duration: 5.3 min) at the KK site. Anatomical scans were acquired using a T1-weighted MPRAGE sequence (GU: TR = 2530 ms, TE = 3.5 ms, FA = 7, voxel size 1 × 1 × 1 mm³; KK: TR = 8.2 ms, TE = 3.7 ms, FA = 8, voxel size 1 × 1 × 1 mm³).

MRI data preprocessing

Data preprocessing was done using functions from FSL (Smith et al., 2004) and AFNI (Cox, 1996); the specific functions are denoted in brackets. Anatomical data was deobliqued (3drefit), oriented into FSL space (RPI) (3dresample) and skull-stripped (3dSkullStrip and 3dcalc). Functional data was also first deobliqued (3drefit) and oriented into FSL space (RPI) (3dresample). The pipeline further consisted of motion correction (MCFLIRT), skull-stripping (3dAutomask and 3dcalc), spatial smoothing (6 mm Gaussian kernel full-width at half-maximum) (fslmaths), grand-mean scaling (fslmaths), registration to the participant's anatomical scan (FLIRT), and normalization to the McConnell Brain Imaging Center NIHPD asymmetrical (natural) pediatric template optimized for ages 7.5–13.5 years (Fonov et al., 2011) (FLIRT), followed by normalization to 2 × 2 × 2 mm MNI152 standard space (FLIRT).

Head-motion and physiological confound mitigation procedure

We used a four-step process to address motion and physiological confounds in the data. First, we used motion estimates derived from the preprocessing in order to exclude participants with excessive head motion; scans were excluded if they exhibited > 4 mm maximum absolute displacement. Second, on the participants who were retained for analysis, we used AROMA, an ICA-based cleaning method (Pruim et al., 2015), which has recently been shown to be most effective in mitigating the impact of head motion (Parkes et al., 2018), and allows for the retention of the remaining ‘true’ neural signal within an affected volume (Kaufmann et al., 2017). AROMA is an automated procedure that uses a small but robust set of theoretically motivated temporal and spatial features (timeseries and power spectrum) to distinguish between ‘real’ neural signals and motion artifacts. We chose a conservative threshold (‘aggressive’) in order to decrease the chance of false positives. In other words, more components are removed as this threshold is more conservative about what is retained. Noise components identified by AROMA were removed from the data. Third, images were de-noised by regressing out the six motion parameters, as well as signal from white matter, cerebral spinal fluid and the global signal, as well as their first-order derivatives (Parkes et al., 2018). While there is currently no gold standard (Murphy and Fox, 2017) regarding the removal of the global signal, we chose to remove it based on recent evidence that it relates strongly to respiratory and other motion-induced signals, which persist through common denoising approaches including ICA and models that attempt to approximate respiratory variance (Power et al., 2018). Motion (defined as each participant’s absolute maximum displacement) was substantially reduced following this procedure (before: 1.28 mm ± 0.85 mm; after: 0.05 mm ± 0.07 mm). As a final step, as described in more detail below, head motion was incorporated into models by removing connections that remained significantly (p < 0.05) associated with z-scored motion before cleaning in a Pearson’s correlation (Rosenberg et al., 2018).

Connectome-predictive modelling

To elucidate how behavioral regulation skills are reflected in children’s whole-brain FC profiles (or ‘connectomes’), and how they vary across children with and without ASD, we utilized a protocol termed Connectome Predictive Modelling (CPM). CPM is an algorithm for building predictive models based on participants’ FC matrices, and for testing these models using cross-validation of novel data. Scripts are written in MATLAB and are freely available at www.nitrc.org/projects/bioimagesuite. The CPM protocol is described in detail in Shen et al., 2017 and has previously been applied to pediatric data sets including data from the ABIDE sample (Lake et al., 2019; Rosenberg et al., 2018, 2016; Finn et al., 2015). We followed the CPM protocol (Shen et al., 2017), as well as recent recommendations for predictive modelling (Scheinost et al., 2019), in calculating each participant’s FC profile, building models of behavioral regulation, and in running the following analyses: (1) a leave-one-out cross-validation to evaluate the potential of models to predict an unseen participant’s score, where N-1 participants are used to build the predictive model and the model is subsequently tested on the left-out participant; (2) a split halves prediction where all available data was randomly split and models built in the first half were used to predict individuals in the second half and vice versa; (3) a site-to-site prediction where models built in the GU site were used to predict individuals in KK site and vice versa. We describe how we calculated FC matrices, built the models and ran these analyses in the following sections.

Calculation of FC profiles

A functionally defined atlas, consisting of 268 cortical and subcortical regions-of-interest (‘nodes’) that cover the whole brain (Shen et al., 2013), was used. For each child, we extracted the timecourse of each node by taking the mean across voxels and a 268 × 268 connectivity matrix was calculated between timecourses of node pairs using Pearson’s correlation followed by Fisher’s Z transformation. Thus, each connection (or ‘edge’) in the matrix represents the strength of FC between two nodes, and the matrix as a whole represents a child’s FC profile or functional connectome.

Building FC models of behavioral regulation

Models were built relating z-scored behavioral regulation subscales (emotional control, shift, inhibit) to FC matrices across participants with and without ASD from both sites. Prior to modeling, effects of motion and site were eliminated from participants’ FC profiles by masking out connections that were significantly (p < 0.05) associated with motion in a Pearson correlation or different between sites in a t-test. 6915 out of a possible 35,778 ( = 268 × 267, adjusted for the diagonal, divided by 2, because matrices are symmetric) nodes were eliminated at this step due to motion, leaving 28,863 valid nodes in the matrix; accounting for site brought this number down to 23100. For model building, each edge in the matrix was correlated with the behavioral regulation measures (again in a Pearson’s correlation), and only significantly correlated edges (p < 0.01) were selected and retained. These selected edges were first separated into positively and negatively associated edges based on the direction of the correlation, as they may be interpreted differently in terms of their functional roles, and then summed for each participant, yielding a single summary FC value per participant for each of the positive and negative edge models. In other words, each participant FC had two summary FC values, one for a network that positively associated with behavioral regulation and one for a network that negatively associated with behavioral regulation. Finally, a predictive model was built that fits a linear regression between each participants’ summary FC value and the behavioral regulation variable of interest (Shen et al., 2017).

Cross-validation: leave-one-out prediction

To evaluate the potential of models to predict an unseen participant’s score, one participant was removed from the dataset and the remaining participants (N-1) were used to build the predictive model. The left-out-participant’s score was predicted based on the N-1 sample’s fit of the linear regression model, and this step was repeated in an iterative manner with a different participant left out in each iteration. Spearman’s rs was used to evaluate model performance i.e. comparing actual to predicted scores because it is less sensitive to the effect of outliers than Pearson’s r and because CPM predictions are best considered relative rather than absolute (Rosenberg et al., 2018, 2016). Only models that showed a significant correlation at p < 0.05 between observed and predicted scores at this step were subjected to follow-up testing.

Evaluation of the predictive model

The predictive potential was assessed by comparison of the predicted and observed scores in the full model, and statistical significance was assessed using permutation testing (5000 iterations). Permutation (i.e., randomization) testing was used to assess significance because the assumption underlying the standard r-to-p conversion employed in the leave-one-out cross-validation (see above) is violated: folds are not independent and thus the number of DOF is over-estimated (Rosenberg et al., 2018, 2016; Scheinost et al., 2019). To perform permutation testing, we randomly shuffled participants’ behavioral scores 5000 times and ran these shuffled values through our prediction pipeline to generate null distributions. P-values associated with each model were based on the corresponding null distribution with the formula p = (1 + the number of permutation rs values greater than or equal to the observed rs value)/5001 (Rosenberg et al., 2018). In other words, the p-value of the permutation test is the proportion of sampled permutations that are greater than the true prediction correlation (Shen et al., 2017).

Cross-validation: split halves prediction

As a further test of the models, data from both sites were randomly split while retaining the same number of TD and ASD participants and the same number of participants from each site. Models were built in the first half and used to predict individuals in the second half and vice versa. Participant characteristics for the two split halves are given in Table 2.

Table 2

Participant characteristics in the two split halves samples. Means and standard deviations (in brackets) are given for the total samples comprised of both TD and ASD participants, which were used to build the models, as well as for TD and ASD participants separately. Motion (mm) refers to the absolute maximum displacement at any timepoint in the resting-state fMRI scan prior to motion mitigation and denoising procedures. n = number of participants; m = male; f = female; L = left-handed; A = ambidextrous; R = right-handed; IQ = Intelligence Quotient; FIQ = full scale IQ; VIQ = verbal IQ; PIQ = performance IQ; SRS = Social Responsiveness Scale. SRS, inhibition, shifting and emotion control are given as T scores. *denotes deviating numbers in the Inhibition models, for which an additional six outliers (>3 SD in score; all ASD participants) were removed. **SRS scores were not available for seven participants in Split Half 1; 1 participant with ASD was removed as an outlier from Split Half 1. † denotes a significant difference between TD and ASD (see Supplementary Table S1 for p-values).

	Split Half 1			Split Half 2
	Total	TD	ASD	Total	TD	ASD
n	138 (135*)	95	42 (39*)	138 (135*)	96	49 (43*)
Age	10.6 (1.4)	10.5 (1.4)	10.6 (1.5)	10.3 (1.3)	10.2 (1.2)	10.3 (1.5)
Sex (m/f)	84/54	53/42	35/7	98/40†	63/33	35/14
Handedness (L/A/R)	8/7/123	5/5/85	3/2/38	13/9/116	8/5/83	5/4/33
Motion	1.2 (0.8)	1.1 (0.8)	1.6 (0.9)	1.3 (0.9)†	1.2 (0.8)	1.8 (0.9)
FIQ	113.6 (12.6)	114.7 (11.3)	108.7 (18.1)	115.1 (14.7)†	117.8 (12.1)	103.7 (15.5)
VIQ	117.6 (14.2)	117.8 (12.9)	112.1 (16.7)	117.8 (14.4)†	120.3 (12.5)	109.4 (16.9)
PIQ	110.4 (12.4)	111 (12.3)	107.7 (16.9)	111.9 (14.3)†	113.5 (12.9)	105.4 (14.3)
SRS**	52.9 (17.6)†	42.7 (5.8)	76.5 (12.1)	53 (15.4)†	44.1 (5.4)	73.1 (10.9)
Inhibition*	48.3 (9.9)†	43.5 (5.4)	59.2 (9.7)	48.9 (10.1)†	44.7 (6.8)	61.7 (8.4)
Shifting	50.4 (14.3)†	43.1 (6.9)	67.8 (12.4)	51.5 (14.1)†	44.3 (7.1)	68.5 (12)
Emotion Control	48.6 (11.4)†	43.1 (6.1)	58.3 (11.99)	48.9 (11.4)†	44.7 (8.3)	59.4 (11.8)

Cross-validation: from site 1 to site 2 prediction

To evaluate the potential of models built in one site to predict an unseen participant’s score from the other site, models were built in the GU site and used to predict individuals in the KK site and vice versa.

Model specificity to behavioral regulation

To evaluate whether our models capture behavioral regulation dimensionally or are driven by the categorical difference in scores due to ASD diagnoses, we examined the relationship between predicted and observed scores for TD children and children with ASD separately. We further examined (a) whether our models could be driven by motion or IQ, (b) how they relate to the core symptoms of ASD, and (c) whether models built for one domain of behavioral regulation were specific to that domain by computing cross-correlations between (i) predicted scores, motion and IQ, (ii) predicted scores and SRS scores, and (iii) predicted and observed scores of different behavioral regulation domains.

Age relationships

The ABIDE dataset only provides T-scores for the BRIEF scales, which are age-normed. For this reason, correlations between age, behavioral regulation and FC models of behavioral regulation were assessed on an exploratory basis, to provide a developmental context.

Results

Sample characteristics

Characteristics for all samples (GU site, KK site, combined sites; split half 1, split half 2) are given in Table 1, Table 2. Scores for all three subscales of behavioral regulation were significantly higher in children with ASD than for TD children in all analyzed samples, reflecting greater challenges with inhibition, shifting, and emotional control (Supplementary Table S1). We further observed significantly higher scores on social responsiveness, reflective of greater ASD symptoms, and several significant differences within and across some of the samples in age, IQ, sex and head motion (Table 2). As expected, the three subscales of behavioral regulation exhibited correlations with each other as well as to social responsiveness and, to a lesser degree, IQ and head motion (Table 3). There were no significant differences in handedness between TD children and children with ASD in any of the samples.

Table 3

Correlations between the three subscales of behavioral regulation, age, sex, motion, IQ and SRS. Results are given for the combined sample and as r-values of bivariate correlations. Motion (in mm) refers to the absolute maximum displacement at any timepoint in the resting-state fMRI scan prior to motion mitigation and denoising procedures. IQ = Intelligence Quotient; FIQ = full scale IQ; VIQ = verbal IQ; PIQ = performance IQ; SRS = Social Responsiveness Scale (total score). *denotes significance at p < 0.05 uncorrected; ** denotes p < 0.0011 (p < 0.05 Bonferroni corrected for 45 comparisons).

	Age	Sex	Motion	FIQ	VIQ	PIQ	SRS	Inhibition	Shifting	Emotion Control
Age		−0.19*	−0.17*	0.05	0.04	0.03	0.03	0.11	0.05	0.01
Sex	−0.19*		−0.07	0.02	−0.05	−0.05	−0.07	−0.01	−0.09	−0.05
Motion	−0.17*	−0.07		−0.18*	−0.10	−0.17*	0.17*	0.14*	0.10	0.11
FIQ	0.05	0.02	−0.18*		0.80**	0.77**	−0.23**	−0.21**	−0.20**	−0.07
VIQ	0.04	−0.05	−0.10	0.80*		0.44**	−0.18*	−0.13*	−0.11	−0.01
PIQ	0.03	−0.05	−0.17*	0.77**	0.44**		−0.18*	−0.14*	−0.16*	−0.05
SRS	0.03	−0.07	0.17*	−0.23**	−0.18*	−0.18*		0.77**	0.83**	0.69**
Inhibition	0.11	−0.01	0.14*	−0.21**	−0.13*	−0.14*	0.77**		0.70**	0.61**
Shifting	0.05	−0.09	0.10	−0.20**	−0.11	−0.16*	0.83**	0.70**		0.79**
Emotion Control	0.01	−0.05	0.11	−0.07	−0.01	−0.05	0.69**	0.61**	0.79**

FC models of behavioral regulation

Significant models were built using negative edges for inhibition (rs = .23, p = 0.0001) and shifting (rs = .19, p = 0.001), using leave-one-out prediction (N-1 at every iteration). Positive edge models were not significant for these measures and neither positive nor negative edge models were significant for emotion control (rs<.05, p > 0.4). The FC model of inhibition was significant by permutation testing (rs = .23, p = 0.037), while the FC model of shifting fell just shy of significance (rs = .19, p = 0.067). As seen in Fig. 1, inhibition revolved around edges in the somato-motor, visual and cerebellar networks and was more posterior, while shifting appeared more focused on edges around the frontoparietal and dorsal attention networks and was more anterior. Both inhibition and shifting included a number of edges connecting with DMN regions as well as the temporal lobe. Note that in the leftmost panels higher rank refers to a lower score, i.e. lower symptoms. In negative edge models, lower FC associates with higher ranked scores.

Fig. 1

CPM models for inhibition (panel A) and shifting (panel B). Models are evaluated using a leave-one-out approach, with a different participant left out in each iteration. The predictive potential is assessed by comparison of the predicted and actual score ranks (left column; inhibition: r = .23; shifting: r = .19) using Spearman’s rank correlation, and statistical significance for the correlation between predicted and observed values is assessed using permutation testing (right column; inhibition: p = 0.037; shifting: p = 0.068). The inhibition model revolved around edges in the somato-motor, visual and cerebellar networks (upper middle column) and was more posterior/inferior, while shifting appeared more focused on edges around the frontoparietal and dorsal attention networks (lower middle column) and was more anterior. Both inhibition and shifting included a number of edges in the default mode network (DMN) and the temporal lobe. The size of the nodes reflects the number of connections the node has to other nodes, with larger nodes being more connected than smaller nodes.

Cross-validation: split halves prediction

Significant models were built for inhibition and shifting using the negative edges in a leave-one-out prediction (N-1 at every iteration) in both split half 1 (inhibition: rs = .26, p = 0.002 and shifting: rs = .32, p = 0.0001) and split half 2 (inhibition: rs = .17, p = 0.049 and shifting: rs = .34, p = 0.00005). The models built in split half 1 further significantly predicted scores in the unseen second half (inhibition: rs = .39, p < 0.000001 and shifting: rs = .19, p = 0.03), and the models built in split half 2 significantly predicted scores in the unseen first half (inhibition (rs = .48, p < 0.000001 and shifting (rs = .19, p = 0.02).

Cross-validation: from site 1 to site 2 prediction

Models could not be built for shifting or inhibition using the negative edges in a leave-one-out prediction (N-1 at every iteration) in the GU site (inhibition: rs = .15, p = 0.17 and shifting: rs = .09, p = 0.4) or in the KK site (inhibition: rs = .02, p = 0.8 and shifting: rs = .03, p = 0.71). Therefore, no cross-prediction from GU to KK and vice versa was attempted.

Model specificity to behavioral regulation

In the combined model, Spearman rank correlations between observed and predicted score ranks for TD children (n = 191) and children with ASD (n = 85) separately were insignificant for the smaller ASD group in both inhibition (rs = .14, p = 0.22) and shifting (rs = .10, p = 0.36), but near significant in the TD group in inhibition (rs = .14, p = 0.053) and significant in shifting (rs = .22, p = 0.002) (Fig. 2). Spearman correlations between predicted shifting or inhibition and IQ or motion, before cleaning or after, were insignificant (rs<.11). Predicted shifting scores associated significantly with total SRS scores (rs = .23, p < 0.001); this association fell below significance in a partial correlation when controlling for diagnosis (rs = .12, p < 0.051). In addition, predicted shifting scores did not significantly associate with observed emotional control (rs = .14, p = 0.018) after correcting for multiple comparisons or inhibition scores (rs = .06, p = 0.33). Similarly, predicted inhibition scores did not associate with observed shifting (rs = .12, p = 0.04) after correcting for multiple comparisons, or emotional control scores (rs = .1, p = 0.09). Associations reported in this section are corrected for the number of comparisons performed (seven per inhibition and shifting, respectively; p < 0.007).

Fig. 2

Spearman rank-correlations between observed and predicted score ranks of inhibition (panel A) and shifting (panel B) for TD children and children with ASD separately, to evaluate whether our models capture behavioral regulation dimensionally or are driven by the categorical difference in scores due to ASD diagnoses. Results were insignificant for the smaller ASD group in both inhibition and shifting, but near significant in the TD group in inhibition and significant in shifting.

Age relationships

Correlations between age and behavioral regulation scores across the combined sample were insignificant (r <.11), as expected, due to T standardization. The FC model for shifting still weakly predicted age (rs=-.12, p = 0.043).

Discussion

Behavioral self-regulation is continually refined across childhood, but many children – especially those with NDCs such as ASD - struggle to regulate their behaviors. In a fully cross-validated, data-driven analysis in a large sample of typically developing children and children with ASD that was compiled across two data collection sites, we identified distributed patterns of FC whose strength predicted individual differences in two behavioral regulation subdomains. These whole-brain network models predicted novel individuals’ inhibition and shifting, but not emotional control scores from resting-state FC data both in a leave-one-out, as well as in a split halves cross-validation, providing evidence that meaningful correlates of behavioral regulation in intrinsic brain patterns exist. Indicating the limitations of this approach, whole-brain network models could, however, not be built within the smaller and less balanced samples collected at the respective sites. We further found that although models captured within-group variation, the model built on shifting also predicted ASD symptoms more generally, although this relationship was no longer significant after taking diagnosis into account. Indicative of the pronounced maturation in FC occurring in relation to behavioral regulation in this age range, the model for shifting also weakly associated with age despite being built on age-normed scores. Overall our results, showing that complex brain network models predict different measures of behavioral regulation across a sample of children with and without ASD, demonstrate that whole-brain FC data can serve as a holistic neural index of inhibition and shifting.

Our findings present a substantial addition to our knowledge on the neural expressions of inhibition and shifting across the spectrum of children with and without ASD. The majority of neuroimaging research on inhibition and shifting has been done in neurotypical adults (Ochsner et al., 2012; Rohr et al., 2016, 2015; Ferri et al., 2016) and adolescents/young adults (Etkin et al., 2010; Kana et al., 2007; Solomon et al., 2009), and little is known about neural alterations underlying shifting and inhibition for children with ASD (Dajani and Uddin, 2015; Yerys et al., 2009) due to varying success of previous research conducted in small sample sizes (Ambrosino et al., 2014; Shafritz et al., 2015; Lee et al., 2009). Encouragingly, our findings are broadly in line with existing literature on shifting and inhibition mechanisms in common ASD comorbidities such as Attention Deficit Hyperactivity Disorder (ADHD), Major Depressive Disorder (MDD) and Obsessive Compulsive Disorder (OCD). This may signify commonalities across these disorders in how inhibition and shifting abilities are reflected in the brain, which bears the implication that if these commonalities were treated in a targeted fashion, inhibition and shifting abilities could be improved across a range of disorders. Both alterations in inhibition and shifting have previously been associated with changes in the DMN - which has hubs in anterior cingulate and ventromedial PFC - in a number of disorders. For instance, activity in regions of the DMN has been shown to be altered in relation to inhibition in ADHD (in a stop signal paradigm (van Rooij et al., 2015)), OCD (during a reward paradigm (Koch et al., 2018)), and MDD (in a cognitive control paradigm (Vanderhasselt et al., 2014)). Similarly, the DMN has been shown to be altered in relation to shifting in Schizophrenia (during reinforcement learning (Waltz et al., 2013)), OCD (Gu et al., 2008) and ADHD (Mulas et al., 2006). In ASD, alterations in FC involving DMN regions have been linked to inhibition (Voorhies et al., 2018) and SRS scores (Jann et al., 2015), which correlate with both inhibition and shifting in our sample.

The divergent neural underpinnings of inhibition and shifting shed further light into children’s functional brain mechanisms. The cerebellum is heavily implicated in ASD (Stoodley, 2016; Becker and Stoodley, 2013; D’Mello et al., 2015) and our group has previously observed that cerebellar FC related to hyperactivity scores in typically developing young children aged 4–7 years old (Rohr et al., 2019), which provide an index of inhibition (Overtoom et al., 2002; Van der Meere et al., 2005). Likewise, somatosensory sensitivities in ASD – like repetitive motor or tactile behaviors – are known to correlate with neural alterations in somato-motor areas (Cascio et al., 2015), and the ability to control repetitive behaviors is linked to inhibition (Schmitt et al., 2018). Individuals with ASD also manifest anomalies in their visual selection and have greater difficulty than neurotypical populations when ignoring specific visual inputs, which relate to alterations in the visual stream (Ronconi et al., 2018). During set-shifting tasks, increased activation in parietal lobes has been reported in individuals with ASD (Schmitz et al., 2006), and our group has previously observed that FC in the dorsal attention network related to attention switching scores in young children (Rohr et al., 2019).

Our work further evidences two major challenges that remain towards achieving one of the primary goals of human neuroimaging - to identify generalizable neuromarkers of clinical utility. First, whole-brain network models could not be built successfully within the smaller and less balanced samples collected at the respective sites. It is a known problem in ASD research that there have been challenges in consistency across the ABIDE sites as differences in demographics and scanner effects in MRI data appear to influence generalizability (King et al., 2019; Chen et al., 2014), meaning that it may be unlikely to build generalizable predictive markers using data from a single site (Scheinost et al., 2019; Nowell et al., 2015). Second, the model built on shifting scores also predicted ASD symptoms generally. While this is to be expected given the known association between ASD core social symptoms and behavioral regulation, it is worthwhile to note that since the relationship between observed and predicted shifting scores remained significant when assessed only within the TD population, they are unlikely to be driven by diagnostic category but rather indicating that the predictions have some specificity to the constructs under investigation. Nonetheless, the correlation between social and behavioral regulation traits makes it more difficult to tease apart the unique and shared features between these constructs.

Individual differences in relation to behavioral regulation have been repeatedly found to be associated with individual features in FC (Rohr et al., 2016; Ferri et al., 2016; Fitzgerald et al., 2019). Taking individual differences into account can help expose the underlying neural substrates of complex cognitive skills, emotions, social competencies and more, and has proven useful in the investigation of both neurotypical (Rohr et al., 2015, 2013; Goldfarb et al., 2016; Vossel et al., 2016) and clinical populations (van Dongen et al., 2015; Nebel et al., 2015; von Rhein et al., 2015). It has been argued that clinical cut-offs for diagnosis may be arbitrary (Uddin et al., 2017; van Dongen et al., 2015), as traits and abilities associated with NDCs such as ASD also exist in the neurotypical population, and fall onto a spectrum or continuum (van Dongen et al., 2015; Matthews et al., 2014). The use of dimensional approaches also allows for more statistical power in studies on neurodevelopmental disorders, which are chronically underpowered due to small sample sizes and challenged by heterogeneity in the populations studied (Uddin et al., 2017; Fair et al., 2012; Nigg, 2005; Sonuga-Barke et al., 2008). Current research suggests that neural patterns associated with abilities that are affected by NDCs are both categorical (i.e. unique to a diagnosis) as well as dimensional (i.e. also present in TD populations) (Elton et al., 2014, 2016). Our findings point to distinct neural mechanisms in the brain subserving the different subtypes of behavioral regulation, which may aid in informing us about options for targeted interventions. They thus highlight possibilities for gleaning insight into how the brain’s functional organization may be associated with cognitive and behavioral issues in children and may serve as a basis for future studies investigating behavioral regulation in neurodevelopmental and other disorders.

Limitations

The current study has several distinct strengths, which include the use of three measures of behavioral regulation in two relatively large, “matched” groups of children with and without an NDC, namely ASD, and novel preprocessing techniques. Both samples were acquired in close spatial and temporal proximity, that is, around the same time and around the same geographic location (Washington, D.C. and Baltimore). The assessment of behavioral regulation is well validated (Gioia et al., 2000) and although parent- and self-reports are subjective, they capture a measure of behavior integrated over a longer time frame than can be observed in a laboratory visit and have better test-retest reliability (Enkavi et al., 2019). At both sites, children were introduced to an MRI simulator first and given the opportunity to familiarize themselves with the experience of undergoing an MRI scan with their eyes open. The denoising methods employed allow for the retention of the remaining ‘true’ neural signal within an affected volume and are in accordance with the latest best practices for reducing motion artifacts (Parkes et al., 2018), which is particularly relevant for a neurodevelopmental study (Pruim et al., 2015; Chen et al., 2015; Zuo and Xing, 2014).

The study also has several weaknesses, including differences within and between the GU and KK sites. These are differences in (i) fMRI data acquisition procedures, (ii) in-/exclusion criteria, (iii) numbers of participants, (iv) TD vs. ASD ratios, (v) sex ratios, (vi) IQs, (vii) BRIEF scores, (viii) comorbidities, and (ix) medications. All of these serve as a reminder of just how divergent ASD populations are, how difficult it is to assess whether a sample reflects the “true patient population” and to what extent research findings are truly generalizable. Considering that these two sites contain a notable number of participants also highlights how common datasets with unbalanced groups are – a known problem when performing cross-validation and assessing prediction performance (Scheinost et al., 2019). Another weakness potentially rests in the removal of motion-associated edges where motion correlates with the behavior that is investigated. While in the combined sample, motion did not significantly associate with the behavioral regulation measures after multiple comparison correction, there remains a possibility that edges that would have associated to, and strengthened, the behavioral regulation models, were removed. The same holds true for site-specific associations. Finally, emotional control scores were not significantly predicted for novel individuals. It is perhaps unsurprising that emotional control results did not follow the same pattern as the other measures of behavioral regulation, as evidence suggests that definitions of emotional control abilities are broad, difficult to dissociate from emotion generation processes, and fairly hard to capture reliably (Gross, 2002; Gross et al., 2011; Mesquita and Frijda, 2011; Preece et al., 2019). Another possibility is that the neural mechanisms of emotional control are simply not reflected in whole-brain FC patterns consistently across individuals.

It should be noted that like most methodological approaches, CPM has both advantages and challenges. One major advantage is that because CPM models are defined and validated with independent data, they promise to improve our ability to uncover generalizable brain-behavior associations (Scheinost et al., 2019; Dubois et al., 2018). A major challenge is that predictive models based on FC will only ever account for a fraction of the variance, because they are limited by how much information the signal can capture as well as the chosen phenotypic measure. However, recent research has highlighted that effect sizes in psychological research are often smaller than previously appreciated, and posited that effect-sizes around those we observe in this study are better understood to indicate “a medium effect that is of some explanatory and practical use even in the short run” (Funder and Ozer, 2019). Predictive models are also bound by linearity assumptions: Linear models built across TD children and children with ASD can capture dimensional associations, but may miss categorical differences which are distinct from dimensional associations (see Elton et al., 2014, 2016 for a discussion on this). Finally, one may consider the current CPM framework as perhaps a bit simplistic in that it yields only one summary value for ‘positive’ and ‘negative’ networks, cannot capture flexible brain network dynamics, and has no ‘blueprint’ for how to tie together predictions on a multitude of behavioral measures.

Future studies could benefit from using data that was obtained while participants perform a task or are under naturalistic viewing conditions such as movie watching. Data obtained while participants perform a task, which adequately captures differences in abilities or skills, has been shown to associate with differences in FC and to lead to better predictive models (Rosenberg et al., 2018; Finn et al., 2017). Showing videos increases young children’s ability to stay still during a scan (Raschle et al., 2009), and may be especially useful for studies in children with NDCs, many of whom evidence attention difficulties in addition to challenges staying still for MRI acquisitions (Rosenberg et al., 2016; von Rhein et al., 2015; Bray et al., 2011a, b, 2013; Keehn et al., 2013; Atkinson and Braddick, 2011). It was further recently shown that individual differences in FC are enhanced during passive viewing, thus facilitating their detection not only through reduced motion but also through the synchronization of hemodynamic fluctuations in large areas of the cortex across participants (Vanderwal et al., 2017). Another improvement could be yielded through the implementation of longer scan imaging times to strengthen the reliability of FC estimates, allow for some data loss in wiggly children and the use of other existing motion mitigation techniques (see e.g. Rohr et al., 2017, 2018; 2019). Finally, harmonized in- and exclusion criteria and scanner and experimental protocols could also aid in providing more comparable FC estimates (Di Martino et al., 2017; Noble et al., 2017).

Conclusions

The characterization of behavioral regulation has been of immense interest to researchers as it matures rapidly in children and is affected in many NDCs and psychiatric conditions with potentially lifelong negative consequences. Yet, it has largely been elusive due to the challenges associated with studying children and the heterogeneity inherent to NDCs such as ASD. In this study, we utilized a data-driven approach to develop objective quantitative FC models that elucidate and predict performance in behavioral regulation subdomains in ASD. We observed both commonalities and differences in the functional organization of inhibition and shifting across TD children and children with ASD, with inhibition relying on more posterior and shifting relying on more anterior brain networks. We also demonstrate the generalizability and trans-diagnostic utility of this approach, as well as its clear limits to date.

Declaration of Competing Interest

None.