This article is mentioned in:

Introduction

Autism spectrum disorder (ASD) is a neurodevelopmental and network-level disorder characterized by impaired development of social interaction, communication skills and a restricted repertoire of activities and interests (1). In 2014, the Centers for Disease Control and Prevention (CDCP) estimated that 1 in 68 children were identified with ASD in the United States, an increase of ~30% compared with previous estimates of 1 in 88 children identified with ASD in 2012 (2,3). In 2014, China's Health and Family Planning Commission reported that an estimated 41,000 children aged 0–6 years were identified with ASD (4). ASD mostly affects language, social communication, motor behaviors and sensory systems (5). Therefore, parents/guardians of children with ASD should consult with pediatric specialists as early as possible, to begin rehabilitation and avoid abuse of diagnosis (6). Currently, the diagnosis of children with ASD is based mostly on behavioral phenotypes (7), which include Autism Diagnostic Interview-Revised (ADI-R) (8), Autism Diagnostic Observation Schedule (ADOS) (9) and the Diagnostic and Statistical Manual of Mental Disorders, Fourth Edition-Text Revision (10). In the United States, only 8% of pediatric specialists capable of conducting ASD screening and diagnosis routinely evaluate toddlers for autism (11). In 2014, China's Health and Family Planning Commission reported that there were ≤100 doctors nationwide trained in diagnosing autism (4). Additionally, the majority of ASD diagnostic tools are not available for use in school-age children. Thus, there is an urgent need for auxiliary diagnosis tools, which may aid to distinguish children with ASD from those typically developing (TD), in particular school-aged children.

Currently, there are several commonly used diagnostic methods and auxiliary diagnostic tools, which include diagnostic methods based on language, behavior, genes, animal models and/or functional magnetic resonance imaging (fMRI) data (12). With the development of interdisciplinary technology, a diagnostic method based on multi-modulus data and machine learning has become popular in recent years (13–16). One of the most effective methods is based on brain frequency and machine learning (14,15,17–19). Zuo et al (20) proposed a concept of brain frequency band and divided the brain frequency band into four frequency sub-bands Slow-2 (0.198–0.25 Hz), Slow-3 (0.073–0.198 Hz), Slow-4 (0.027–0.073 Hz) and Slow-5 (0.01–0.027 Hz). Other studies also demonstrated that the method of diagnosing autism based on brain frequency was possible (21–24). These studies also demonstrated that different brain frequency bands revealed different properties and physiological functions and they hypothesized that these could be the ‘fingerprints’ of neuronal activity (25). However, these investigations based on brain frequencies demonstrated limited use for the classification of children with ASD and those TD. In addition, the accuracy of classification decreased when the fMRI data was obtained from a single institution (19). Chen et al (19) used fMRI data from multiple institutions to distinguish ASD from TD and demonstrated a 74% diagnostic accuracy when using two frequency bands compared with a 69% diagnostic accuracy when using only one frequency band. The heterogeneity and complexity of ASD has impaired the research of ASD in children (26). Therefore, methods based on one frequency band and traditional machine learning cannot be used for the auxiliary diagnosis of ASD. Full brain frequency bands contain more information, however they also have some noise, such as physiological noise (20). The stacked auto-encoders (SAE) deep learning framework is a neural network consisting of multiple layers and it can extract important features from input data (27). In the current study, an advanced deep learning method based on SAE, full frequency bands and multi-institutional fMRI data was used to distinguish school-aged children with ASD from those TD. The current study examined large fMRI datasets of 198 school-aged children from three different institutions obtained from the Autism Brain Imaging Data Exchange (ABIDE) II database (fcon_1000.projects.nitrc.org/indi/abide/index.html) (28) and demonstrated relatively high classification accuracy compared with previous studies (16,19,29–31).

Materials and methods

Data collection

Original fMRI datasets of 198 school-aged children (5–12 years) were downloaded from the ABIDE II database, which was released in June 2016 and enables unrestricted use for non-commercial purposes. The current study was approved by the Ethics Committee of Nanchang University (Nanchang, China).

Patient charts

The current study used original fMRI datasets of school-aged children obtained from the New York University Langone Medical Center (New York, NY, USA), Georgetown University (Washington, DC, USA) and the Kennedy Krieger Institute (Baltimore, MD, USA). For each participant, a high-resolution structural MRI sequence was acquired using different MRI scanners with different parameters, as well as a resting-state fMRI sequence using different MRI scanners with different parameters. The names of the three different institutions and the different fMRI scanning parameters are summarized in Table I.

Table I. fMRI scanning parameters and settings.

Table I.

fMRI scanning parameters and settings.

Institution	MRI manufacturer	TR (msec)	TE (msec)	Voxel (mm3)	Slices	Scan time (min)	Eyes
Georgetown University	Siemens AG, Munich, Germany	2,000	30	3.00×3.00×2.50	43	5.23	Open
Kennedy Krieger Institute	Philips Medical Systems, Inc., Bothell, WA, USA	2,500	30	2.67×2.67×3.00	47	6.67	Open
New York University	Siemens AG, Munich, Germany	2,000	30	3.00×3.00×3.00	34	6.00	Open
Langone Medical Center

[i] fMRI, functional magnetic resonance imaging; MRI, magnetic resonance imaging; TR, repetition time; TE, echo time.

All participants were diagnosed according to the score obtained in both the ADI-R and the ADOS. The TD participants were screened during clinical interviews conducted by child psychiatrists and details of the diagnostic procedures and questionnaires are listed in the ABIDE II database (28). According to age, sex and full-scale intelligence quotient (FIQ), 117 participants with ASD (minimum age, 5.22 years; maximum age, 11.99 years; mean age [standard deviation (SD)], 9.32 [5.62] years; FIQ [SD], 107.87 [15.26]) and 81 age- and FIQ-matched TD participants [minimum age, 5.88 years; maximum age, 11.95 years; mean age (SD), 9.67 (5.31); mean FIQ (SD), 109.40 (15.12)] were selected from the datasets obtained from New York University Langone Medical Center, Georgetown University and the Kennedy Krieger Institute. Demographic data for the selected participants are summarized in Table II.

Table II. Demographic data of ASD and TD participants.

Table II.

Demographic data of ASD and TD participants.

Participant group	Patients (n)	Male/female ratio	Minimum age, years	Maximum age, years	Mean age, years (SD)	Mean FIQ (SD)
ASD	117	93/24	5.22	11.99	9.32 (5.62)	107.87 (15.26)
TD	81	57/24	5.88	11.95	9.67 (5.31)	109.40 (15.12)

[i] ASD, autism spectrum disorder; TD, typically developing; FIQ, full-scale intelligence quotient; SD, standard deviation.

Resting-state fMRI data preprocessing

Resting-state fMRI preprocessing was performed using the multivariate exploratory linear optimized decomposition into independent components (MELODIC; version 3.14), a tool of FMRIB Software Library (FSL) software packages (32). Data preprocessing included head motion correction slice timing correction, removal of non-brain tissues, and spatial normalization into MNI152 space using nonlinear registration and spatial smoothing using 6 mm full width half maximum Gaussian kernel. Currently, the prediction of autism using brain networks involves examining functional connectivity in 10 functional networks using independent component analysis (33). A previous study identified 30 independent components (ICs) corresponding to 10 previously described functional networks (34). By contrast, several studies revealed that the number of ICs increases to 50, additionally the time of fMRI data decomposition increases exponentially and the dimension of obtaining brain frequency signals also increases (26,35). High-dimensional brain signals will directly make the neural node of the deep learning algorithm increase, which will also increase the training time of the system. Therefore, each participant's data were decomposed into 30 ICs using MELODIC (36).

Frequency selection and normalization

Following IC decomposition using MELODIC version 3.14, the time serial, the power spectrum and the spatial map are generated for each IC (Figs. 1–3). MELODIC can generate the ICs on one spatial map followed by the relevant time serial of the IC and the power spectrum of the time serial. The power spectrum is an expression of brain frequency power distribution. The full brain frequency was divided into 5 different bands: Slow-6 (0–0.01 Hz), Slow-5 (0.01–0.027 Hz), Slow-4 (0.027–0.073 Hz), Slow-3 (0.073–0.198 Hz) and Slow-2 (>0.198 Hz). The power spectra of 30 ICs for each participant contain the information of each sub-frequency band in the network. In order to get more information, all ICs were selected in a data-driven manner although some ICs may contain useless information such as noise. All ICs with full brain frequency were selected according to the frequency ranges, although some ICs may contain low frequency drifting and noise.

Figure 1. Time series decomposition following IC analysis. Red line represents the time course of IC11. IC, independent component.

Figure 3. Spatial map of an independent component. Spatial map of IC11. IC, independent component.

For each participant, a time-series of all ICs were reshaped into a feature array with the same dimension. The feature array was normalized in the range of (0,1) using Max-Min normalization before it was input into the SAE-based classifier.

SAE-based classification

In general, a sparse auto-encoder (AE) is used for feature extraction and dimensionality reduction. The sparse AE is a framework of unsupervised learning method, and SAE is a neural network consisting of multiple layers of sparse AE, where the output of each layer is connected to the input of each successive layer (35). A softmax classifier is used for phase classification, a supervised learning approach (37). SAE and softmax classifier were integrated into the proposed computer-aided diagnosis method, as shown in Fig. 4. AE is a key substructure of the SAE, with a symmetrical architecture (Fig. 5). The outputs of each layer in SAE were wired to the inputs of the successive layer. As shown in Fig. 6, the structure of the enhanced SAE proposed in the current study contains 10 layers: An input layer, 8 hidden layers and an output layer. The 8 hidden layers were divided into 4 groups, each with 2 hidden layers. The numbers of nodes in each of the 4 groups were as follows: 1,280, 320, 80 and 20.

Figure 4. Proposed classification system. fMRI, functional magnetic resonance imaging; MELODIC, multivariate exploratory linear optimized decomposition into independent components; ICA, independent component analysis; SEA, stacked auto-encoders; TD, typically developing; ASD, autism spectrum disorder.

Figure 5. Structure of a sparse auto-encoder.

Figure 6. Schematic architecture of an enhanced stacked auto-encoder.

The i-th activation unit a(l)i of the l-th layer can be represented by the sl−1 units of (l−1)-th layer through the weight parameter set w, a bias term set b and an activation function f as shown below in Equations 1 and 2:

[n=number of units in the (l-1)-th layer, f(•) was chosen as the sigmoid function].

Given a training dataset with m samples, {(x(1), y(1)),…,(x(m), y(m))} where y(i) was the category label of i-th sample, y(i) ϵ {1, 2,…,k}, for the first layer, the activation units were the input data (i.e. a(1)i=xi). AE enabled the output values to gradually approach the input values by applying a back propagation algorithm. Therefore, the goal of AE was to minimize the distance between the inputs and outputs to learn w and b. So the cost function was computed as shown below in Equation 3:

[m was the number of training samples, nl (nl=3) was the number of layers, sl was the number of units in the l-th layer and λ was the weight decay parameter]. In equation 3, the first term was the mean square error between input and output values, representing the quality of learning, while the second term was a regularization term that tends to decrease the magnitude of the weights and to prevent the learning from overfitting.

To ensure that the outputs of the hidden layer are desirably sparse, a sparsity constraint must be introduced into the cost function to control the learning process. So the average activation of the j-th hidden unit was defined as (with respect to input x) as shown below in Equation 4:

If ρ̂j was small (close to 0.01), most of the units of the hidden layer should be inactive. To force ρ̂j to be equal to a very small value ρ, the sparsity penalty term was designed based on the concept of Kullback-Leibler divergence, as shown below in Equation 5:

Subsequently, the sparsity penalty term was incorporated into the cost function, as shown below in Equation 6:

(β was the weight of sparsity penalty term). During the process of learning, the cost function Jsparse (W, b) was minimized by updating w and b. As presented in (31), limited-memory BFGS (L-BFGS) was a suitable optimization algorithm for updating w and b in the process of back propagation. In the current study, L-BFGS was adopted, and the softmax regression model was used to classify school-aged children with ASD and TD. After data preprocessing, the time serial matrix was input into the multi-layer SAE, and the output values from the SAE were input into the softmax classifier. The current study used MELODIC (version 3.14), an analysis tool used to decompose each participant's fMRI data into ICs. The time-series was extracted from each IC, and used to generate a vector time-series of each participant with 30 time-series. Similarly, the system was used to obtain 198 vector time-series, and all vector time-series were normalized to the same dimension in order to build a matrix. The matrix was divided into training matrix and test matrix according to the different fold cross-validation (CV). The training matrix was input into the SAE for optimization training. Finally, the trained multi-layer SAE and softmax classifier were used to classify participants of the testing matrix. In this cycle, multiple test results can be obtained. In the current study, the classification accuracy of all subjects was computed using a CV. For example, in the 198-fold CV, the original samples are randomly partitioned into 198 groups, with each group including one participant. Of the 198 groups, 197 groups (197 participants) were chosen as the training set, and the remaining one participant was used as the test data. The validation was repeated 198 times such that each participant was used exactly once as the test data. Averaging the 198 results during the validation produced the final estimation. The classification statistical indicators included true positive (TP), false negative (FN), true negative (TN) and false positive (FP) values. The parameters of accuracy, sensitivity, specificity, positive prediction value (PPV) and negative prediction value (NPV) were computed using the following formulae: Accuracy=(TP+TN)/(TP+TN+FP+FN) ×100% (7); sensitivity=TP/(TP+FN) ×100% (8); specificity=TN/(TN+FP) ×100% (9); PPV=TP/(TP+FP) ×100% (10); and NPV=TN/(TN+FN) ×100% (11).

The sensitivity measured the proportion of positive results that were correctly identified. Specificity measured the proportion of negative results that were correctly identified, or the percentage of TD children who were correctly identified as TD. The whole process was performed using Dell Precision T5810 (CPU: Xeon (R), 3.5 GHz; MM: 128 GB; OS: 64-bit Windows 10 professional edition) and Matlab R2012a (7.14.0.739, 64 bits) was used to run the multi-layer SAE and softmax regression algorithm. The specific algorithm settings including, the weight decay, sparsity penalty term, weight of sparsity term and the number of nodes, are summarized in Table III.

Table III. Algorithm parameter values used in the 8 hidden layers.

Table III.

Algorithm parameter values used in the 8 hidden layers.

Hidden layer	Weight decay	Sparsity penalty term	Weight of sparsity penalty term	Number of node
No. 1, 2	1×10−8	0.01	3	1,280
No. 3, 4	1×10−8	0.01	3	320
No. 5, 6	1×10−8	0.01	3	80
No. 7, 8	1×10−8	0.01	3	20

Results

In the current study, the classification accuracy was evaluated using CV. The classifier assigned a label (ASD=1, TD=2) to each participant, and the total accuracy was computed for each category. The accuracy, specificity, sensitivity, PPV and NPV for each CV are summarized in Table IV.

Table IV. Cross-validation classification accuracy analysis.

Table IV.

Cross-validation classification accuracy analysis.

Parameter	11-fold (%)	33-fold (%)	66-fold (%)	99-fold (%)	198-fold (%)	Average, % (SD)
Accuracy	94.19	93.69	98.99a	97.47	96.97	96.26 (2.02)
Specificity	87.04	87.23	100a	96.30	97.53	93.62 (5.43)
Sensitivity	99.15a	97.86	98.29	98.29	96.58	98.03 (0.84)
PPV	91.70	91.97	100a	97.46	98.26	95.88 (3.40)
NPV	98.60a	96.60	97.59	97.50	95.18	97.09 (1.15)

a Maximum value for each parameter analyzed. PPV, positive prediction value; NPV, negative prediction value; SD, standard deviation.

With the same 198 subject data and 198-fold CV, the SAE, probabilistic neural network (PNN) and support vector machine (SVM) algorithms achieved 96.97, 84.58 and 83.69% accuracy, respectively (Table V). Therefore, the proposed SAE with 8 hidden layers of the current study was more accurate than other machine learning algorithms. These results demonstrated that the SAE algorithm was a more effective method for the analysis and processing of data in the current study. In addition, the 66-fold CV results revealed that the proposed SAE algorithm produced a better balance between the classification accuracy and computation time compared with SAE algorithms with different numbers of hidden layers (data not shown). The average accuracy obtained using the SAE algorithm with 6, 8 and 10 hidden layers was 93.69, 98.99 and 95.90%, respectively (data not shown). However, the computation time taken for SAE algorithms with 6, 8 and 10 numbers of hidden layers was 70,018, 68,588 and 75,891 sec, respectively (data not shown). In conclusion, increasing the number of hidden layers did not improve the accuracy of the SAE algorithm, however the running time increased. These results suggest that the SAE algorithm with 8 hidden layers can balance the classification accuracy and running time.

Table V. Comparison of classification accuracy of three different algorithms.

Table V.

Comparison of classification accuracy of three different algorithms.

	Stacked auto-encoders (%)	Probabilistic neural network (%)	Support vector machine (%)
Accuracy	96.97	84.58	83.69

The current study used a novel deep learning method with full brain frequency, to generate a higher diagnostic accuracy (96.26%), with increased sensitivity (98.03%) and specificity (93.62%). Furthermore, the accuracy classification analysis of the current study was compared with previously published studies (Table VI). The results demonstrated that the classification accuracy was ~15% higher in the current study compared with previous studies, where only one or two frequency bands were used (19,33,38–40).

Table VI. Comparison between the classification accuracy analysis obtained in the current study and previously published studies.

Table VI.

Comparison between the classification accuracy analysis obtained in the current study and previously published studies.

Author, year	Functional networks	Accuracy (%)	Sensitivity (%)	Specificity (%)	PPV (%)	NPV (%)	(Refs.)
Current study		96.26	98.03	93.62	95.88	97.09	–
Uddin et al, 2013	Salience	78.00	75.00	80.00	79.00	76.00	(33)
	Central executive	58.00	55.00	60.00	58.00	57.00	(33)
	Posterior default mode	63.00	65.00	60.00	62.00	63.00	(33)
	Ventral default mode	60.00	55.00	65.00	61.00	59.00	(33)
	Anterior default mode	63.00	50.00	75.00	67.00	60.00	(33)
	Motor	68.00	60.00	75.00	71.00	65.00	(33)
	Visual association	65.00	65.00	65.00	65.00	65.00	(33)
	Dorsal attention	73.00	75.00	70.00	71.00	74.00	(33)
	Primary visual	73.00	60.00	85.00	80.00	68.00	(33)
	Frontotemporal	68.00	60.00	75.00	71.00	65.00	(33)
Anderson et al, 2011 (subject <20 years)		89.00	83.00	95.00	–	–	(38)
Chen et al, 2016 (NYU dataset)		80.00	95.00	72.00	–	–	(19)
Chen et al, 2016 (UM dataset)		84.00	69.00	94.00	–	–	(19)
Chen et al, 2016 (USM dataset)		80.00	80.00	80.00	–	–	(19)
Iidaka, 2015 (50-fold CV)		90.30	92.30	88.40	88.30	92.40	(40)
Iidaka, 2015 (10-fold CV)		86.90	85.90	87.80	87.00	86.70	(40)
Iidaka, 2015 (2-fold CV)		77.20	75.00	79.30	77.50	76.90	(40)

[i] PPV, positive prediction value; NPV, negative prediction value; Refs., references; NYU, New York University; UM, University of Michigan; USM, Universiti Sains Malaysia; CV, cross-validation.

Discussion

In the current study, the proposed SAE algorithm-based classification system was generated using all frequency bands using multi-institutional fMRI datasets to produce high classification accuracy to distinguish school-aged children (5–12 years) with ASD from those TD. The present study examined an auxiliary method to distinguish ASD with high sensitivity and specificity, which was considered as an entirely hypothesis-free and data-driven method. In a previous study, ASD (12–42 years) and TD (8–39 years) participants were classified with relatively low accuracy (38). As ASD is an early onset disorder with variable developmental trajectory, Uddin et al (33) used a network-based classification method to predict children (7–13 years) with ASD and those TD. Although Uddin et al (33) examined 10 whole-brain functional networks (salience, central executive, posterior default mode, ventral default mode, anterior default mode, dorsal attention, motor, visual association, primary visual and frontotemporal), the salience network was the only network, which demonstrated good classification accuracy. In the current study, several other brain networks were integrated, and the classification accuracy was increased compared with single brain network classification as demonstrated by Uddin et al (33). Another study classifying ASD and TD participants through a default network only achieved a relatively low accuracy compared with the current study (39). Although previous studies gained some success in mechanism research, previous studies used relatively low numbers of participants and the results were based on single brain networks only. However, the data-driven method used in the current study used a bigger dataset, and attained a higher accuracy by not only integrating single brain networks, but also measuring the joint effect of several brain networks. In 2015, Iidaka (40) used a PNN approach on 640 participants to distinguish between ASD and TD, with relatively high accuracy.

Following a comparison with previously published studies using PNN and SVM methods, the current study revealed that the SAE algorithm achieved a higher classification accuracy compared with traditional machine learning algorithm. In 2016, Chen et al (19) used SVM to distinguish between ASD and TD, and obtained relatively low classification accuracy. Although Chen et al (19) and Iidaka (40) used larger datasets, which were obtained from the same database, the 2 classification algorithms (SVM and PNN) demonstrated less accuracy compared with the current study. The differences in methodology between the SAE and SVM/PNN algorithms may account for the differences observed in classification accuracy of the above-mentioned studies and the current study.

In conclusion, a deep learning based approach was used to successfully distinguish 198 school-aged children with ASD from those TD, using a data-driven method with high accuracy of 96.30%. Taken together, these results demonstrate the potential clinical application of ASD diagnostic tools in school-aged children with ASD.

Acknowledgements

The authors would like to thank the researchers and funding agencies that contributed to the ABIDE II database.

Funding

The current study was funded by grants from The Natural Science Foundation of China (grant nos. 61662047 and 61463035).

Availability of data and materials

All datasets used and/or analyzed during the current study are available from the corresponding author in reasonable request. In addition, the original functional magnetic resonance imaging datasets can be obtained from the Autism Brain Imaging Data Exchange II database (fcon_1000.projects.nitrc.org/indi/abide/index.html).

Authors' contributions

ZX designed the experiment, analyzed the functional magnetic resonance imaging (fMRI) data, interpreted the results and prepared the manuscript. JW and CW contributed significantly to the conception and design of the work, and were involved in the collection and interpretation of data. XY prepared all the figures in the manuscript and detected the stacked auto-encoders algorithm. NJ performed the preprocessing of the fMRI data. All authors read and approved the final manuscript.

Ethical approval and consent to participate

The current study was approved by the Ethics Committee of Nanchang University. In addition, the ABIDE II database allows unrestricted use for non-commercial purposes.

Patient consent for publication

Not applicable.

Competing interests

The authors declare that they have no competing interests.

References