Abstract
Background:
To develop a new functional magnetic resonance image (fMRI) network inference method, BrainNET, that utilizes an efficient machine learning algorithm to quantify contributions of various regions of interests (ROIs) in the brain to a specific ROI.
Methods:
BrainNET is based on extremely randomized trees to estimate network topology from fMRI data and modified to generate an adjacency matrix representing brain network topology, without reliance on arbitrary thresholds. Open-source simulated fMRI data of 50 subjects in 28 different simulations under various confounding conditions with known ground truth were used to validate the method. Performance was compared with correlation and partial correlation (PC). The real-world performance was then evaluated in a publicly available attention-deficit/hyperactivity disorder (ADHD) data set, including 134 typically developing children (mean age: 12.03, males: 83), 75 ADHD inattentive (mean age: 11.46, males: 56), and 93 ADHD combined (mean age: 11.86, males: 77) subjects. Network topologies in ADHD were inferred using BrainNET, correlation, and PC. Graph metrics were extracted to determine differences between the ADHD groups.
Results:
BrainNET demonstrated excellent performance across all simulations and varying confounders in identifying the true presence of connections. In the ADHD data set, BrainNET was able to identify significant changes (p < 0.05) in graph metrics between groups. No significant changes in graph metrics between ADHD groups were identified using correlation and PC.
Conclusion:
We describe BrainNET, a new network inference method to estimate fMRI connectivity that was adapted from gene regulatory methods. BrainNET out-performed Pearson correlation and PC in fMRI simulation data and real-world ADHD data. BrainNET can be used independently or combined with other existing methods as a useful tool to understand network changes and to determine the true network topology of the brain under various conditions and disease states.
Impact statement
Developed a new functional magnetic resonance image (fMRI) network inference method named as BrainNET using machine learning.
BrainNET out-performed Pearson correlation and partial correlation in fMRI simulation data and real-world attention-deficit/hyperactivity disorder data.
BrainNET does not need to be pretrained and can be applied to infer fMRI network topology independently on individual subjects and for varying number of nodes.
Introduction
The brain is a complex interconnected network that balances segregation and specialization of function with strong integration between regions, resulting in complex and precisely coordinated dynamics across multiple spatiotemporal scales (Sporns, 2018). Connectomics and graph theory offer powerful tools for mapping, tracking, and predicting patterns of disease in brain disorders through modeling brain function as complex networks (Fornito et al., 2015). Studying brain network organization provides insight in understanding global network connectivity abnormalities in neurological and psychiatric disorders (Avena-Koenigsberger et al., 2018). Several studies suggest that pathology accumulates in highly connected hub areas of the brain (Buckner et al., 2009; Crossley et al., 2014) and that cognitive sequelae are closely related to the connection topology of the affected regions (Warren et al., 2014). An understanding of network topology may allow prediction of expected levels of impairment, determination of recovery following an insult, and selection of individually tailored interventions for maximizing therapeutic success (Fornito et al., 2016). A large number of network inference methods are being used to model brain network topology with varying degrees of validation. A recent study (Smith et al., 2011) evaluated some of the most common methods, including correlation, partial correlation (PC), and Bayes NET, to infer network topology using simulated resting-state functional magnetic resonance image (rs-fMRI) data with known ground truth and found that performance can vary widely under different conditions.
Development of statistical techniques for valid inferences on disease-specific group differences in brain network topology is an active area of research (Kim et al., 2015, 2019). Machine learning methods have been used in neuroimaging for disease diagnosis and anatomic segmentation (Murugesan et al., 2018; Zaharchuk et al., 2018). Brain Network Construction and Classification (BrainNetClass) and GraphVar toolboxes provide a full pipeline from network construction to classification. BrainNetClass comprises various fMRI network inference methods such as correlation, PC, and higher order functional connectivity for brain network inference followed by feature extraction for machine learning model development and testing (Waller et al., 2018; Zhou et al., 2020). Very few studies have attempted to apply machine learning methods on direct time series of fMRI to infer brain networks (O'Neill et al., 2017; Pellegrini et al., 2018; Williams and Henson, 2018; Zaharchuk et al., 2018). Recent work in machine learning approaches for inference of gene regulatory networks (GRN) has demonstrated excellent performance (Camacho et al., 2018; Finkle et al., 2018; Turki et al., 2016). Interestingly, these same approaches to GRNs can be used to infer brain networks. In this study, we describe a new network inference method called BrainNET, inspired by machine learning methods used to infer GRN (Irrthum et al., 2010).
Yan et al. (2018) devised a bidirectional long short-term memory (LSTM) deep learning network (Full-BiLSTM) to effectively learn the periodic fMRI brain status changes using both past and future information for each brief time segment. They then fused them to form the final output by taking a dynamic functional connectivity matrix calculated using the sliding window approach as input (Yan et al., 2018). Higher order functional connectivity was developed by Chen et al. (2016) by taking dynamic relationships between the brain regions to infer network topology. Yu et al. (2017) proposed a novel method using connectivity weighted sparse representation to construct optimal brain functional networks from rs-fMRI data. The method has taken advantage of both Pearson's correlation and sparse representations, which are the two most commonly used brain network modeling approaches. This ensures the construction of more biologically meaningful brain networks by a unified framework that integrates connectivity strength, group structure, and sparsity. Yu et al. (2017) used l1-norm regularized linear regression or sparse representation, which learns a linear relationship, while BrainNET considers the nonlinear relationships. The abovementioned methods, including BrainNetClass and GraphVar, focus on using machine learning methods for computer-aided diagnosis by predicting the cognitive metrics or classifying a group of the subjects, whereas BrainNET uses machine learning for inferring the networks directly.
Validation of BrainNET was performed using fMRI simulations with known ground, as well as in real-world attention-deficit/hyperactivity disorder (ADHD) fMRI data sets. In this study, publicly available resting-state fMRI simulated data (Smith et al., 2011) were used to validate BrainNET's ability to infer networks. The real-world performance of BrainNET was then evaluated in a publicly available data set of ADHD. ADHD is one of the most common neurodevelopmental disorders in children with significant socioeconomic and psychological effects (Hilger and Fiebach, 2019; Lin et al., 2014). It can be difficult to diagnose due to the overlapping nature of symptoms, with resultant diagnostic errors and overprescribing of medications due to misdiagnosis (Saeed, 2018). ADHD has widespread but often subtle alterations in multiple brain regions affecting brain function (Cortese et al., 2012; Sidlauskaite et al., 2015). Neuro Bureau, a collaborative neuroscience forum, has released fully processed open-source fMRI data “ADHD-200 preprocessed” from several sites (Bellec et al., 2017; Milham et al., 2012) providing an ideal data set to test the BrainNET model and compare its performance with standard correlation and PC, the most widely used methodology to infer brain networks using fMRI data.
Materials and Methods
Data sets
MRI simulation data
Open-source rs-fMRI simulation data representing brain dynamics were used to validate the BrainNET model (Smith et al., 2011). The data were simulated based upon the dynamic causal modeling fMRI forward model, which uses the nonlinear balloon model for vascular dynamics, in combination with a neural network model (Smith et al., 2011). The open-source data set has 28 simulations, each including simulated data for 50 subjects with a varying number of nodes and several confounders (e.g., shared input between the nodes, varying fMRI session lengths, noise, cyclic connections, and hemodynamic lag variability changes). Additional details on the simulations can be found in the original study (Supplementary Table S1) (Smith et al., 2011).
ADHD data
Preprocessed rs-fMRI data were obtained from the ADHD-200 database (
In our study, subjects identified with “no naive medication” status, or questionable quality on rs-fMRI data was excluded. The remaining subjects were age-matched between the groups resulting in 135 typically developing children (TDC) (mean age: 12.00, males: 83), 75 ADHD inattentive (ADHD-I) (mean age: 11.46, males: 56), and 93 ADHD combined (ADHD-C) (mean age: 11.86, males: 77) subjects. Mean time series from 116 regions of interests (ROIs) in the automated anatomic atlas (AAL) atlas (Tzourio-Mazoyer et al., 2002) were extracted using the NILEARN package (Abraham et al., 2014).
BrainNET model development
The objective of BrainNET is to infer the connectivity from fMRI data as a network with N different nodes in the brain (i.e., ROIs), where edges between the nodes represent the true functional connectivity between nodes. At each node, there are measurements from m time points
Our method assumes that fMRI measurement of BOLD (blood oxygen level dependent) activation at each node is a function of each of the other nodes' activation with additional random noise.
For the jth node with m time points, a vector can be defined denoting all nodes except the jth node as
where ɛj is random noise specific to each node j . We further assume that function ƒj () only exploits the data of nodes in x−j that are connected to node j . The function ƒ j () can be solved in various ways in the context of machine learning. Since the nature of the relationship between different ROIs in the brain is unknown and expected to be nonlinear (Stam and Reijneveld, 2007), we choose a tree-based ensemble method as it works well with a large number of features with nonlinear relationships and is computationally efficient. We utilized extremely randomized trees (ERT), an ensemble algorithm similar to random forest, which aggregates several weak learners to form a robust model. ERT uses a random subset of predictors to select divergences in a tree node and then selects the “best split” from this limited number of choices (Praagman, 1985). Finally, outputs from individual trees are averaged to obtain the best overall model (Petralia et al., 2015). BrainNET infers a network with N different nodes by dividing the problem into N different subproblems, and solving the function ƒ j () for each node independently as illustrated in Figure 1. The steps are listed below:

Schematic overview of the BrainNET model. For N nodes in fMRI data (X), each node will have m time points such that
. Each node's time series (xn
) is predicted from all other nodes time series (x-n
) using the ERT regressor. Node importance of each node for predicting the target node is extracted and populated in the importance matrix. The average of the upper and lower triangle of the matrix is thresholded at (1/Num of Nodes) to obtain an adjacency matrix representing the network topology. ERT, extremely randomized trees; fMRI, functional magnetic resonance images. Color images are available online.
For j = 1 to N nodes Fit the ERT regressor with all the node data, except the jth node, to find the function fj
that minimizes the following mean squared error:
Extract the weight of each node to predict node j,
where wn
is the weight of node to predict node j and n = 1 to N. Append the weight values to the importance matrix.
The importance score for each node (Node j ) to predict (Node i ) is defined as the total decrease in impurity due to splitting the samples based on Node j (Praagman, 1985). The GINI index is used here as the measure of impurity. Let “S” denote a node split in the tree ensemble and let (SL, SR) denote its left and right children nodes. Then, the decrease in impurity ΔImpurity(S) from node split “S” based on Node j to predict Node i is defined as follows:
where SL and SR are left and right splits and NP, NL, NR are number of samples reaching parent, left, and right nodes, respectively. Let
be the number of ensembles, which uses ROI
j
for splitting trees. Then, the importance score for Node
j
for predicting Node
i
is calculated as the average of node impurities across all trees, that is, importance of ROI
ji
where T is the number of trees in the ensemble.
Importance values extracted using a typical random forest model can be biased in the presence of two or more correlated features since the model will randomly assign importance to any one of the equally important features without any preference (Strobl et al., 2007). This problem is avoided by using the ERT regressor. The comparison between ERT, random forest, and baseline LASSO in inferring network topology with the simulation data is provided in Supplementary Figure S1 and Supplementary Table S2.
The importance of each node to predict all other node time series is extracted from the model and an N × N (where N is the number of nodes) importance matrix is generated with the diagonal equal to zero. Each row of the importance matrix represents normalized weights of each node in predicting the target node. The extracted adjacency matrix is affected in two ways. First, due to the row-wise normalization, the upper triangular values of the importance matrix are not same as the lower triangle values. We therefore take the average of the upper triangle and the lower triangle of the matrix to make it symmetric to determine the presence of connection between the nodes. This procedure does not allow directionality of the connections to be determined. The comparison between ERT, random forest, and baseline LASSO in inferring network topology with the lower and upper triangle, averaged (symmetrized) and nonaveraged (nonsymmetrized), is provided in Supplementary Figure S1 Second, again because of the row-wise normalization, the sum of each row in the importance matrix is one. Since the importance values are normalized with respect to number of nodes in the analysis, we used a threshold that is inversely proportional to the number of nodes (i.e., threshold = 1/number of nodes) in the network to produce a final adjacency matrix representing the network topology. The selection of threshold is not based on statistical theory and it is not made to keep the false positive rate below a nominal level, but it results in a dynamically changing threshold based on the number of nodes in the network.
Analysis
Evaluation of inference methods on simulation data
Evaluation of inference methods on simulation data using c-sensitivity
The network topology was inferred using BrainNET, correlation, and PC. The network topology inferred by correlation and PC method may vary drastically based on the values used to threshold connectivity matrix. Hence, we evaluated the ability of the inference methods based on BrainNET, correlation, and PC to detect the presence of connection between the nodes in terms of c-sensitivity. C-sensitivity quantifies how well the true positives (TPs) are separated from the false positives (FPs) by measuring the fraction of TPs that are estimated with a higher connection strength than the 95th percentile of the FP distribution. C-sensitivity is a measure of success in separating true connections from FP connections and it is calculated by counting the number of TPs above the 95th percentile of FPs and then divided by the total number of TPs (Smith et al., 2011).
Effects of simulation parameters such as TR (repetition time), number of nodes, noise, hemodynamic response function standard deviation, shared inputs, bad ROIs, backward, strong, and cyclic connections, and strong inputs on c-sensitivity of the inference methods were evaluated using the mixed-effects model with a random effect for the simulation to control for the effect of the specific generating model. The mixed-effects models were fit across subjects under different simulations to analyze the effects of simulation parameters on c-sensitivity of the inference methods. The parameter estimates from each regression were then summarized across subjects in terms of their effect size (Smith et al., 2011).
Evaluation of inference methods on simulation data using threshold
Thresholding can be applied to suppress spurious connections that may arise from measurement noise and imperfect connectome reconstruction techniques and to potentially improve statistical power and interpretability (Fornito et al., 2016). However, based on the threshold value, the connection density of each network inferred by correlation and PC may vary from network to network after the threshold has been applied. Using a less stringent lower threshold value results in higher FP values (lower sensitivity) and a more stringent threshold results in higher false negatives (lower specificity). This can lead to wide variability in computed graph metrics, as they are typically susceptible to the number of edges in a graph. Identifying an appropriate threshold to infer the underlying brain network topology is critical.
Hence, we evaluated the specificity, sensitivity, and accuracy of correlation and PC under varying thresholds. The results show that the network topology inferred using correlation and PC method may vary drastically based on the threshold values (Fig. 2). An optimum threshold for correlation methods can be very difficult to find in real-life experimental data. However, given the ground truth for the simulation data, we calculated the optimum threshold values for the correlation methods and compared their performance at optimum threshold with BrainNET. The optimal thresholds are defined where the correlation and PC methods performed with the best sensitivity and accuracy. It is important to note that BrainNET is not optimized on these simulation data, and the threshold is based on the number of nodes inferred in the network. Specificity, sensitivity, and accuracy of correlation and PC at threshold values of 30% (Corr30, PC30) and optimum (Corropt, PCopt) values are estimated and compared with BrainNET. We further evaluated the specificity, sensitivity, and accuracy of the Corropt and PCopt with BrainNET for each simulation. Similar to the c-sensitivity, we calculated the effect of simulation parameters on the sensitivity and specificity of BrainNET, Corropt, and PCopt.

Sensitivity analysis for correlation and PC. Average sensitivity (true-positive rate) and specificity across 28 simulations for the correlation and PC methods are plotted as a function of threshold ranging between 0% and 100%. Optimum threshold is found using simulation ground truth at 20% and 16% for correlation and PC, respectively. PC, partial correlation. Color images are available online.
Evaluation of inference methods on ADHD data
BrainNET, correlation, and PC were applied to the real-world ADHD data to evaluate whole-brain network changes in ADHD subtypes (i.e., ADHD-C, ADHD-I) compared with TDC. Mean time series from 116 ROIs in the AAL atlas (Tzourio-Mazoyer et al., 2002) were extracted using the NILEARN package (Abraham et al., 2014). The BrainNET model was applied to extract an importance matrix for each subject. The importance matrix was then thresholded at 1/number of nodes (e.g., 1/116 for the AAL atlas regions) to obtain an adjacency matrix for each subject (BN). Functional network connectivity was calculated between the 116 ROIs using correlation and PC. The connectivity matrices are thresholded at a threshold of 20% and 30% (Corr20, Corr30, PC20, and PC30) (no optimum threshold for real-world experimental data). Graph theoretic metrics were extracted using each of these methods for each group. Network differences between the three groups TDC, ADHD-I, and ADHD-A were then computed using t-tests on the graph metrics. Site effects and effects of age and handedness were removed using the ComBat multisite harmonization method (Yamashita et al., 2019), an effective harmonization technique that both removes unwanted variation associated with the site and preserves biological associations in the data (Fortin et al., 2018).
Graph metrics
Graph theoretical metrics representing global and local characteristics of network topology were used to compare between the groups in the ADHD data set. The GRETNA MATLAB toolbox (v2.0,
Node metrics
The nodal shortest path length (NSPL) is defined as the shortest mean distance from a particular node to all other nodes in the graph. Shorter NSPL represents greater integration (Sporns, 2018). The BC measures a node's influence in information flow between all other nodes (Wang et al., 2010). BC quantifies the influence of a node and is defined as the number of shortest paths passing through it.
Global metrics
Network efficiency is a more biologically relevant measure representing the ability of the network to transmit information globally and locally. Networks with high efficiency, both globally and locally, are said to be economic small-world networks (Achard and Bullmore, 2007). The density of the graph is defined as the ratio of number of connections in the network to the number of possible connections in the network. Average clustering is the fraction of a node's neighbors that are neighbors of each other. The clustering coefficient of a graph is the average clustering coefficient over all nodes in the network. Networks with high clustering coefficient are considered locally efficient networks. CPL is the average shortest path length between nodes in the graph, with a minimum number of edges that must be traversed to get from one node to another. CPL indicates how easily information can be transferred across the network (Sporns, 2018).
Experimental Results
Simulation data
Evaluation of inference methods using c-sensitivity
BrainNET performed significantly better than correlation (p < 0.001) and equivalent to PC (p > 0.05) methods with c-sensitivity of 79.53%, 59.82%, and 75.75%, respectively, across 28 simulations (Fig. 3A and Supplementary Table S3).

C-sensitivity. Boxplots of c-sensitivity for BrainNET, correlation, and PC (left). The effects of different simulation parameters on the c-sensitivity, sensitivity, and specificity of inference methods using mixed-effects model. The color bar represents the effect size of each simulation parameter on the c-sensitivity, sensitivity, and specificity of inference methods (right). Color images are available online.
The study on the effect of simulation parameters on the c-sensitivity of inference methods showed that increasing the number of nodes and session duration does not have much effect on any of the inference methods (Fig. 3B). Having shared inputs between the nodes affected the correlation and PC method more drastically than the BrainNET method. Selection of bad ROIs with mixed time series between them affected all the inference methods negatively, however, selection of bad ROIs with randomly random time series mixed between them did not affect the inference methods drastically. The presence of cyclic and backward connection between the nodes affected the correlation and PC methods but not BrainNET. The presence of only one strong input affected the performance of the BrainNET method but not the other methods. In summary, BrainNET performance was robust under various confounding factors but prone to selection of inaccurate ROIs with mixed time series between them and networks with only one strong input. Both PC and correlation methods were affected by shared inputs between the nodes, selection of inaccurate ROIs, and backward, cyclic, and stationary/nonstationary connections between the nodes.
Evaluation of inference methods using threshold
The accuracy, sensitivity, and specificity for each method across all 28 simulated data sets were estimated at thresholds of 30% and at optimum value for correlation and PC, and at a threshold of 1/number of nodes for BrainNET (Fig. 4). BrainNET achieved higher accuracy and specificity at a threshold of 30% compared with the Corr30 and PC30 method, as shown in Figure 4. PCopt achieved a slightly higher accuracy than BrainNET across 28 simulations, but no significant difference in terms of specificity and accuracy even at its optimum (p > 0.05) (Fig. 5). As expected, the specificity and sensitivity of correlation and PC methods vary with threshold and it will be difficult to find an optimum threshold in a real-life data set. BrainNET showed more robust performance with little variance across the simulation compared with other methods.

Evaluation of inference methods under varying thresholds. Boxplots of accuracy (left), sensitivity (middle), and specificity (right) across 28 simulations for correlation and PC for optimum and 30% threshold (Corropt, Corr30, PCopt and PC30), and BrainNET. * Represents statistically significant differences from BrainNET performance. Color images are available online.

Comparison of correlation (Corropt) and PC (PCopt) at their optimum threshold to BrainNET. Accuracy (left), sensitivity (middle), and specificity (right) for correlation, BrainNET, and PC for 28 simulations. Sensitivity, specificity, and accuracy are all robust across different simulation cases, while PC and correlation methods show fluctuations even with their optimal threshold for functional connectivity. Color images are available online.
The study on effect of simulation parameters on sensitivity and specificity after threshold for BrainNET, PCopt, and Corropt showed that the sensitivity of the inference methods does not get affected by the simulation parameters (Fig. 3B). The specificity of the correlation and PC methods was negatively affected even at their optimum threshold in the presence of backward connections. It is important to note that the specificity and sensitivity were calculated after we thresholded the connectivity matrices from each of these methods (1/number of nodes for BrainNET and optimum threshold for others). After the threshold, BrainNET's performance became robust and consistent across all the simulation parameters.
Evaluation of inference methods on ADHD data
BrainNET was able to identify significant changes (p < 0.05) in global network efficiency, network density, CPL, BC, and shortest path in the ADHD data. Correlation and PC were not able to detect significant changes in any of the whole-brain analyses (Corr20, Corr30, PC20, and PC30).
TDC and ADHD
Statistical analysis of the BrainNET adjacency matrix demonstrated a significant decrease in global network efficiency, an increase in CPL, and an increase in the shortest path length in the right medial temporal gyrus in ADHD compared with TDC (Fig. 7A). While the analysis of the correlation adjacency matrix did not show any significant changes, the PC30 demonstrated a trending increase in CPL in ADHD compared with TDC (p = 0.07). BC and node level local efficiency did not show any changes between the groups in any of the three methods.
TDC and ADHD-I
Statistical analysis of the BrainNET adjacency matrix demonstrated a significant decrease in global network efficiency, a decrease in density, an increase in CPL, and an increase in the shortest path length in the right superior orbital, right Heschl's gyrus, and right medial temporal gyrus nodes in the ADHD-I group compared with TDC (Fig. 7A). The correlation method did not show any relationship between the groups. No relationship was found in other graph metrics studied in any of the methods.
TDC and ADHD-C
Statistical analysis of the BrainNET adjacency matrix demonstrated a significant decrease in density. No significant relationship was found in any other graph metrics for any of the three methods.
ADHD-I and ADHD-C
Statistical analysis of the BrainNET adjacency matrix demonstrated a significant decrease in global network efficiency, a decrease in density, an increase in CPL, and an increase in shortest path length of the right olfactory node (Fig. 7A). A significant increase in BC of the right precuneus node in the ADHD-I group compared with ADHD-C was observed for BrainNET (Fig. 7B). No relationship was found in other graph metrics studied for any of the methods.
Discussion
BrainNET was developed to infer brain network topology using ERT (Geurts et al., 2006). The ERT regressor is used to develop a tree-based ensemble model to predict each node's time series from all other node time series. The tree-based ensemble methods are ideal for inferring complex functional brain networks as they are efficient in learning nonlinear patterns even where there are a large number of features (Wehenkel et al., 2017). The importance matrix is then thresholded to generate an adjacency matrix representing the fMRI topology. The BrainNET model is applicable to both resting-state and task-based fMRI network analysis. It can be easily adapted to data sets with varying session lengths and can be used with different parcellation schemes. A unique feature of the BrainNET approach is that it is implemented at the subject level. It does not need to be trained on big data sets as it infers the network topology based on each individual subject's data.
BrainNET inference of network topology in simulated fMRI data
Evaluation of inference methods using c-sensitivity
BrainNET demonstrated excellent performance across all the simulations and varying confounders. It achieved a significantly higher c-sensitivity than correlation (p < 0.05) and equivalent to PC (p = 0.38) (Fig. 2A). BrainNET performance remained high in the simulations across varying session lengths, number of nodes, neural lags, cyclic connections, and changing number of connections. BrainNET performed weakest in simulations with one primary strong external source around the network. This causes every node to be highly correlated with other nodes and it becomes very difficult to distinguish direct from indirect connections (Smith et al., 2011). It is important to highlight that this kind of one strong external input just for one node is highly unlikely in real-life scenarios. BrainNET, similar to PC and correlation methods, was affected by selection of bad ROIs with time series mixed between them. In this simulation, there are 10 nodes, and each node shares a relatively small amount of the other node time series in a proportion of 0.8:0.2. Since the features have shared data between the nodes in this simulation, it limits discrimination of true connectivity between nodes. The leakage of data between nodes can be minimized in fMRI analysis by selecting independent regions using functionally derived parcellation or methods such as independent component analysis.
One concern with this approach is that as the number of nodes increases, the threshold [1/(number of nodes)] similarly decreases and may result in increased FPs at this low threshold value. The study on effect of number of nodes on the c-sensitivity of BrainNET shows that c-sensitivity of BrainNET does not get affected by the number of nodes (Fig. 2B). This can be interpreted that the ability of BrainNET to distinguish between TP and FPs increases with the increasing number of nodes and the corresponding lower threshold values do not necessarily affect its inference.
Shared inputs between the nodes can be thought of as distinct sensory inputs that feed into one or more nodes. These shared inputs between the nodes could be deleterious if not modeled (Smith et al., 2011). BrainNET is robust to the shared inputs between the nodes, whereas c-sensitivity of PC and correlation is negatively affected (Fig. 2B). The performance on varying connection strengths over time was tested by simulations of stationary/nonstationary connection strengths between the nodes. BrainNET was least affected by nonstationary connection strengths between the nodes (Fig. 2B). The robust performance of BrainNET in simulations with increasing number of nodes, TR, shared inputs, and backward, cyclic, and nonstationary connections represents a promising aspect of the BrainNET method for inferring brain network topology in real-life experimental data (Fig. 2B).
Evaluation of inference methods using thresholds
In this study, we compared the performance of correlation and PC in inferring underlying network topology at optimum threshold values estimated using ground truth (Corropt and PCopt). We also compared the performance of these methods against BrainNET at 30% threshold (Corr30 and PC30) (Fig. 4). BrainNET performed significantly better than PC30, Corr30, and Corropt (p < 0.05). At its optimum threshold, PCopt performed relatively equivalent to BrainNET in terms of accuracy. However, PCopt showed decreased specificity with increasing number of nodes (sim1–4) and the presence of nonstationary and backward connections (sim13 and sim22) (Fig. 5). Nonstationary connections represent the varying strengths of connections between nodes, which are believed to be similar to those at the neuronal level and being studied in fMRI. The higher sensitivity and lower specificity of PCopt represent higher numbers of FP connections, which will affect the statistical power of group analysis. The results show that the performance of the PC and correlation method varies under different thresholds and that BrainNET had a better performance than these methods even in their optimum (PCopt and Corropt) (Fig. 4). The study on the effect of the number of nodes on c-sensitivity of BrainNET shows robust performance across all the confounders.
A major strength of the BrainNET approach is that it provides a unique threshold to determine the true network topology. In correlation-based approaches, there is no defined correlation cutoff to determine the true network topology. Instead, multiple approaches are used, or multiple thresholds are applied to generate different networks. Typically, the network cost has been used to define the cutoff value for defining true connections in correlation-based approaches (Supekar et al., 2009). Multiple costs are then applied to generate multiple instances of the network topology, and analyses are performed to determine the variation in network metrics across these costs, or variation in group differences across thresholds (Achard and Bullmore, 2007). The BrainNET approach provides a single threshold obviating the need for these imprecise and convoluted thresholding approaches.
Evaluation of inference methods on ADHD data
Global metrics
Previous studies have shown that ADHD is often associated with changes in the functional organization of the brain and lower network efficiencies in ADHD (Lin et al., 2014; Sidlauskaite et al., 2015). BrainNET was effective in identifying the subtle changes in the ADHD subjects and supports the notion that the functional organization of brain changes in ADHD (Fig. 6) by identifying statistically significant changes in graph metrics between ADHD subjects and TDC.

Global graph metrics. The probability density functions and boxplots of global graph metrics with significant changes (p < 0.05) between the groups, ADHD (both ADHD-I and ADHD-C), ADHD-I, ADHD-C, and TDC. ADHD, attention-deficit/hyperactivity disorder; ADHD-C, ADHD combined; ADHD-I, ADHD inattentive; CPL, characteristic path length; TDC, typically developing children. For ease of reading and color images, the figure can be viewed online.
Our results demonstrate that there is a decrease in density and network efficiency, and an increase in CPL in ADHD compared with TDC. A decrease in density suggests that the number of connections is decreased in ADHD compared with TDC. This can be interpreted as an increase in the cost of wiring in the brain. The increase in CPL and the decrease in network efficiency are expected, given that there is a decrease in density suggesting that there is increased difficulty in transferring information across the brain in ADHD. The observed abnormalities in global network topology were identified in ADHD-I, but not in participants with ADHD-C compared with TDC, however, changes between the ADHD-C and ADHD-I were observed. The differential changes observed between the ADHD subtypes may reflect clinical distinctions between the inattentive and combined subtypes of ADHD. Further investigations may shed light on detailed brain-behavior phenotype associations in this neuropsychiatric disorder (Barber et al., 2015; Qian et al., 2019).
Local metrics
BrainNET identified increased NSPL in ADHD-I compared with ADHD-C, suggesting lesser integration of the prefrontal cortex (PFC) in ADHD-I. The PFC is a part of the default mode network (DMN) and plays a crucial role in regulating attention, behavior, and emotion, with the right hemisphere specialized for behavioral inhibition (Arnsten, 2009). The DMN refers to the brain circuitry that includes the medial PFC, posterior cingulate, precuneus, and the medial, lateral, and inferior parietal cortices (Weyandt et al., 2013). These results support previous studies demonstrating that ADHD is associated with structural changes and decreased function of the PFC circuits, especially in the right hemisphere (Arnsten, 2009). BrainNET also demonstrated that BC of the right precuneus, also a part of DMN, was increased in ADHD-I compared with the ADHD-C group. This suggests increased influence of the precuneus in ADHD-I (Fig. 7B). Abnormalities within the DMN have also been found in children in previous studies with ADHD and especially changes in centrality of the right precuneus, which is an important discriminatory feature for classifying ADHD-I and ADHD-C (dos Santos Siqueira et al., 2014).

Node level graph metrics. Nodes with significant increases in NSPL in ADHD-I compared with ADHD-C (orange) and in ADHD-I compared with ADHD-C (red) are plotted in the left. Nodes with significant increases in betweenness centrality in ADHD-I compared with ADHD-C are plotted on the right. NSPL, nodal shortest path length. Color images are available online.
Our results also show that the NSPL of the right Heschl's gyrus and right medial temporal gyrus is increased in the ADHD-I group compared with TDC. The NSPL of the olfactory cortex was increased in ADHD-I compared with ADHD-C (Fig. 7A). Deficits in olfactory function are found in neurodegenerative and neuropsychiatric disorders and represent a topic of interest in ADHD (Ghanizadeh et al., 2012). Increased NSPL was found in the right olfactory region in ADHD-I compared with ADHD-C, suggesting lesser integration. Deficits in olfactory ability have been linked to impulsive tendencies within the healthy population and have discriminatory features in identifying people at risk of impulse-control-related problems, supporting the planning of early clinical interventions (Herman et al., 2018). Further studies are required to investigate whether the functional network topology can be used as a biological marker for early diagnosis, treatment, and prognosis of ADHD.
It is important to note that the proposed method measures nonlinear relationships, while the correlation methods measure linear relationships, which may have resulted in the lower performance of correlation in inferring nonlinear brain dynamics. Although PC performed relatively similarly to BrainNET in the simulation data, it did not achieve statistical significance in the ADHD data. This may be due to the FPs identified reducing the statistical power of the analysis. BrainNET can be added to the standard inference methods such as PC and correlation methods, by using a mask derived from the BrainNET importance matrix and applying to the correlation matrix. The output from this combined method will have nodes determined by BrainNET, with Pearson correlation values assigned between the connections. This will avoid using arbitrary thresholds, increase the specificity of the standard inference methods by adding nonlinearity and allowing analysis of connectivity changes between nodes, which cannot be performed with an adjacency matrix derived only from BrainNET.
Limitations
BrainNET takes relatively longer to infer the adjacency matrix than the correlation method. BrainNET took ∼3 sec per subject, whereas the correlation and PC methods just took 0.001 and 9.3 sec, respectively. Longer running time makes BrainNET challenging to apply for voxel-wise analysis.
Conclusion
We describe BrainNET, a new network inference method, to estimate fMRI connectivity that was adapted from gene regulatory methods. We validated the proposed model on ground truth simulation data (Smith et al., 2011). BrainNET outperformed Pearson correlation in terms of accuracy and sensitivity across simulations and various confounders such as the presence of cyclic connections, and even with truncated fMRI sessions of only 2.5 min. We evaluated the performance of BrainNET on the open-source “ADHD 200 preprocessed” data from Neuro Bureau. BrainNET was able to identify significant changes in global graph metrics between ADHD groups and TDC, whereas correlation and PC were unable to find any differences. BrainNET can be used independently or combined with other existing methods as a useful tool to understand network changes and to determine the true network topology of the brain under various conditions and disease states.
Footnotes
Author Contributions
Gowtham Krishnan Murugesan conceptualized and designed the work, analyzed and interpreted the data, and wrote the article. Dr. Joseph Maldjian and Dr. Won Hwa Kim provided expert knowledge and mentorship to develop the method. Ben Wagner contributed in developing fMRI analysis. Chandan Ganesh, Sahil Nalawade, and Dr. Elizabeth Davenport contributed in reviewing the article.
Author Disclosure Statement
No competing financial interests exist.
Funding Information
No funding was received for this article.
Supplementary Material
Supplementary Figure S1
Supplementary Table S1
Supplementary Table S2
Supplementary Table S3
References
Supplementary Material
Please find the following supplemental material available below.
For Open Access articles published under a Creative Commons License, all supplemental material carries the same license as the article it is associated with.
For non-Open Access articles published, all supplemental material carries a non-exclusive license, and permission requests for re-use of supplemental material or any part of supplemental material shall be sent directly to the copyright owner as specified in the copyright notice associated with the article.
