Discriminative analysis of non-linear brain connectivity in schizophrenia: an fMRI Study

Background: Dysfunctional integration of distributed brain networks is believed to be the cause of schizophrenia, and resting-state functional connectivity analyses of schizophrenia have attracted considerable attention in recent years. Unfortunately, existing functional connectivity analyses of schizophrenia have been mostly limited to linear associations. Objective: The objective of the present study is to evaluate the discriminative power of non-linear functional connectivity and identify its changes in schizophrenia. Method: A novel measure utilizing the extended maximal information coefficient was introduced to construct non-linear functional connectivity. In conjunction with multivariate pattern analysis, the new functional connectivity successfully discriminated schizophrenic patients from healthy controls with relative higher accuracy rate than the linear measure. Result: We found that the strength of the identified non-linear functional connections involved in the classification increased in patients with schizophrenia, which was opposed to its linear counterpart. Further functional network analysis revealed that the changes of the non-linear and linear connectivity have similar but not completely the same spatial distribution in human brain. Conclusion: The classification results suggest that the non-linear functional connectivity provided useful discriminative power in diagnosis of schizophrenia, and the inverse but similar spatial distributed changes between the non-linear and linear measure may indicate the underlying compensatory mechanism and the complex neuronal synchronization underlying the symptom of schizophrenia.


INTRODUCTION
Schizophrenia, which is characterized by delusions, auditory hallucinations, and impairments in memory, attention, and executive function, is one of the most devastating, cryptic, and costly psychiatric disorders (van Os et al., 2010;Sui et al., 2012). It brings not only great suffering to patients but also significant costs to society. Traditionally, the diagnosis of schizophrenia depends on the observation of psychiatric symptoms and longitudinal courses, while modern diagnoses of such psychiatric disorders require objective neurological measures (Kendler, 2009;Insel et al., 2010;Shen et al., 2010). The underlying etiology and mechanisms of schizophrenia are still unclear, but it is believed that the dysfunctional integration of distributed brain networks leads to this mental disease (Friston and Frith, 1995;Andreasen et al., 1998). Thus, based on fMRI data, functional connectivity research on schizophrenia has attracted considerable attention in recent years (Camchong et al., 2011;Pettersson-Yeo et al., 2011;Bassett et al., 2012;Fornito et al., 2012;Liu et al., 2012;Tu et al., 2012;Vertes et al., 2012;Zalesky et al., 2012). Furthermore, the introduction of multivariate pattern classification to behavioral and cognitive neuroscience increases the potential of functional connectivity in clinical diagnoses of this mental disease (Shen et al., 2010;Fan et al., 2011).
Functional connectivity, defined as statistical associations between remote neurophysiological events, aims to characterize the communication between different brain regions (Friston et al., 1993). Most functional connectivity analyses use temporal correlations or covariance to examine the simultaneous coupling between two time series. Fluctuations in the blood oxygen level dependence (BOLD) signal have attracted attention since the 1990's (Biswal et al., 1995(Biswal et al., , 1996. In contrast to functional connectivity, which need previous acknowledge about the disease and did not estimate the potential functional connectivity of other ROIs (Lowe et al., 1998(Lowe et al., , 2000Xiong et al., 1999;Cordes et al., 2000;Hampson et al., 2002), whole brain functional connectivity analyses of schizophrenia which require less field specific knowledge have attracted considerable attention in recent years Liu et al., 2006;Lynall et al., 2010;Alexander-Bloch et al., 2013).
Previous whole brain functional connectivity analyses of schizophrenia which primarily used temporal correlation or covariance as the measure limited to consideration of only the linear associations Liu et al., 2006;Lynall et al., 2010;Alexander-Bloch et al., 2013). Because correlation calculations in the full-lag space are computationally expensive (Cecchi et al., 2007), correlation related studies have prevalently used the zeroth lag correlation. The hemodynamic response can significantly reduce the computational expense for its limited duration (Li et al., 2009), but presents another problem, namely, the hemodynamic response function (HRF) varies across subjects and brain regions (Buckner et al., 1998;Miezin et al., 2000;Lee et al., 2001;Saad et al., 2001), and the zeroth lag correlation is sensitive to changes in regional HRFs. Specifically, the simple correlation (zeroth lag) is sensitive to the lag between time series that cannot sufficiently depict the functional interactions of the human brain (Smith et al., 2011). Thus, we believe that the underlying neural activity cannot be accurately reflected by the zeroth lag correlation-or covariance-based methods.
Investigating functional connectivity in schizophrenia from the frequency domain is an alternative of the temporal correlation (Fallani et al., 2010;Salvador et al., 2010). Although the spectral analogs of functional connectivity such as coherence fixed some of the problems (i.e., lag between time series) faced by the simple correlation methods, coherence explores only linear relationships between time series (Sun et al., 2004;Smith et al., 2011). However, the non-linearity of the HRF has been reported by several studies (Buxton et al., 1998;Friston et al., 2003;de Zwart et al., 2009;Daunizeau et al., 2012). As an output of HRF, the BOLD signal also has non-linear properties (Vazquez and Noll, 1998;Xie et al., 2008a,b;Zhang et al., 2008). Notably, non-linear relationships between time series extracted from resting state BOLD signals have also been confirmed (Lahaye et al., 2003).
The first goal of this study is to evaluate the discriminative power of non-linear functional connectivity in schizophrenia which may have potential application in diagnosis of neuropsychiatric disorders. The second goal is to investigate the changes of the non-linear functional connectivity in schizophrenia which may contribute to the pathophysiology of schizophrenia. Aiming at exploring non-linear associations in schizophrenia, we present a novel method that introduces the extended maximal information coefficient (eMIC) to construct functional connectivity of the whole brain.
The eMIC represents the difference between the MIC which is a statistical method for detecting various associations between pairs of variables in large data sets (Reshef et al., 2011) and the square of the Pearson correlation coefficient (PCC). In a previous study, the square of PCC was proved to be equal to the Hilbert-Schmidt Independence Criterion (HSIC) and gave excellent performance in feature selection (Song et al., 2007). By applying the MIC, PCC and eMIC with fMRI data, we investigated the discriminative power of non-linear and linear functional connectivity. Based on the classification result, we further evaluated the changes and spatial distribution of linear and non-linear functional connectivity in schizophrenia.

Participants
The fMRI data used in this study were acquired from 64 participants who were all right-handed native Chinese speakers. Participants consisted of 32 patients suffering from schizophrenia and 32 healthy controls. All the patients were recruited from outpatient departments and inpatient units at the Department of Psychiatry, Second Xiangya Hospital of Central South University in Changsha, China, between March 2006 andOctober 2007, andsatisfied the Diagnostic andStatistical Manual of Mental Disorders, Fourth Edition (DSM-IV;American Psychiatric Association, 1994) criteria. Five patients were medication-free, while the others accepted antipsychotic drugs at the time of image acquisition. The healthy subjects were recruited by advertisements, and matched to the patients on age and gender. None of them had major head trauma, history of alcohol or drug dependence, or history of neurological disorder. Written informed consents were obtained from all the subjects who took part in this study. This study was approved by the Medical Research Ethics Committee of the Second Xiangya Hospital, Central South University. Details about the participants were displayed in Table 1.

Imaging protocol
In the experiments, subjects were instructed simply to keep their eyes closed, to relax, remain awake and perform no specific cognitive exercise. After each session, subjects were asked whether they were awake in the previous session and all the subjects confirmed. MRI scans were performed with a 1.5T GE Signa System (GE Signa, Milwaukee, Wisconsin, USA). To reduce the head movements, subjects' heads were fixed by using foam pads with a standard birdcage head coil. The functional MRI images were collected by using a gradient-echo echo planar imaging sequence. The imaging parameters were as follows: repetition time/echo time = 2000/40 ms, thickness/gap = 5/1 mm, field of view = 240 × 240 cm, flip angle = 90 • , matrix = 64 × 64, and slices = 20. Each functional resting state session lasted ∼6 min, and 180 volumes were obtained.

Data preprocessing
Data preprocessing was performed by using statistical parametric mapping software package (SPM2, Wellcome Department of Cognitive Neurology, Institute of Neurology, London, UK, http://www.fil.ion.ucl.ac.uk/spm/). In each subject, the first 5 volumes of the scanned data were discarded for magnetic saturation effects. The remaining 175 volumes were corrected by registering and reslicing for head motion. Next, the volumes were normalized to the standard echo planar imaging template in the Montreal Neurological Institute (MNI) space. The resulting images were spatially smoothed with a Gaussian filter of 8 mm full-width half-maximum kernel, detrended to abandon linear trend and then temporally filtered with a Chebyshev band-pass filter (0.01-0.08 Hz). The registered fMRI volumes with the MNI template were further divided into 116 regions according to the automatic anatomical labeling atlas (Schmahmann et al., 1999;Tzourio-Mazoyer et al., 2002).The atlas divides the cerebrum into 90 regions (45 in each hemisphere) and divides the cerebellum into 26 regions (nine in each cerebellar hemisphere and eight in the vermis). All region of interest masks were generated using the free software WFU_PickAtlas (version 3.0, http://www.ansir. wfubmc.edu) (Maldjian et al., 2003). Regional mean time series were obtained for each individual by averaging the functional MRI time series over all voxels in each of the 116 regions. Each regional mean time series was further corrected for the effects of head movement by regression on the time series of translations and rotations of the head estimated in the course of initial movement correction by image realignment. The residuals of these regressions constituted the set of regional mean time series used for functional connectivity analysis (Achard et al., 2006).
We evaluated functional connectivity between each pair of regions using the MIC, PCC, and eMIC. Thus, for each subject, we obtained three resting state functional networks captured by three 116 × 116 symmetric matrixes respectively. Removing 116 diagonal elements, we extracted the upper triangle elements of the functional connectivity matrix as classification features, i.e., the feature space for classification was spanned by the (116 × 115)/2 = 6670 dimensional feature vectors.

MIC and emic
In this section, we provide a brief description of the MIC and eMIC for detecting the association between two time series. Two time series can be viewed as a set of ordered data pairs. The MIC between a set of ordered pairs will not change if the rank of pairs is disturbed but the relative ranks of the x-and y-values do not change. If two variables are related to each other, then a grid can be drawn on the scatterplot of the two variables that encapsulates that relationship. Based on this concept, this method investigated all the grids up to the maximal grid resolution, which depends on the sample size (Reshef et al., 2011). The MI is defined as follows: where H(X) and H(Y) are the marginal entropies, H(Y |X ) and H(Y |X ) are the conditional entropies, and H(X, Y) is the joint entropy of X and Y.
Formally, let G x × y be all the possible partitions with rows x and y columns (width of rows are different; width of the columns are different, too) of the scatterplot for the ordered pairs of two vectors. As the partitions were not equal, there are many possible partitions with x rows and columns of the scatterplot, let I g denote the MI for one possible partition with x × y grids that are applied to the ordered samples of the two vectors. For fair comparison between grids of different resolution, the MI values of different partitions with x × y grids of scatterplot for the ordered pairs of two vectors are normalized to the interval [0, 1]. The m x × y is defined as MIC is the maximal value of m x × y over the ordered pairs (x, y), x ≤ n, y ≤ n, n is the length of the vectors (i.e. number of factors in the vector). In practice, to accelerate computation, it is not necessary to compute m x × y over all (x, y), x = n, y = n. Alternatively, we can compute the m x × y over all (x, y) xy < B.
Empirically, B is defined as B = n 0.6 . We can then define Additionally, the eMIC can be defined as Refer to (Song et al., 2007;Reshef et al., 2011) for more details of the MIC and HSIC.

Identification of Features with High Discriminative Power
Due to the noise, registration error, and inter-individual anatomical differences, only a small number of the 6670 features are highly discriminative. To achieve good discriminative performance as well as resistance to noise and individual disparity, the first step of constructing the classification model was selecting those features with high discriminative power to construct the feature space for classification. The discriminative power of a feature can be quantitatively measured by its relevance to classification (Guyon and Elisseeff, 2003). Here, we used the Kendall tau rank correlation coefficient (Kendall and Jean, 1990) which provides a distribution free test of independence between two variables. Let x 1 , y 1 , x 2 , y 2 , . . . , x N , y N be a set of samples of the joint random variables X and Y respectively. Any pair of samples x i , y i and x j , y j are said to be concordant if the ranks for both elements agree. The Kendall tau coefficient is defined as: of this sample (+1 for controls and −1 for patients). The Kendall tau correlation coefficient of the functional connectivity feature i can be defined as Because the samples in the same groups are neither concordant pairs nor discordant pairs, the relationship between two samples that belong to the same group is not considered, number of total simple pairs The n c and n d are the numbers of concordant and discordant pairs between the two groups, respectively. For a pair of two-observation data sets {x ij , y j } and {x ik , y k }, it is neither concordant pair nor discordant pair if y i = y k ; it is a concordant pair when where sgn( ) is a signum function. Correspondingly, it is a discordant pair when Thus, a positive correlation coefficient τ i represents that the ith functional connectivity exhibits significant decrease in the patient group compared to the control group, while a negative τ i represents that the ith functional connectivity increases in the patient group. The discriminative power was defined as the absolute value of the Kendall tau correlation coefficient. Then we ranked all the features according to their discriminative power and selected those with correlation coefficient over a threshold as the final feature set for classification.
Since a leave-one-out cross-validation strategy was introduced to estimate the generalization ability of the classifiers (see below) and the training data set for feature ranking is slightly different in each iteration of the cross-validation, the first selected features differed slightly from iteration to iteration. Therefore, the contribution of different regions to classification was not evenly distributed, and some regions formed many highly discriminating functional connections with other regions, while some did not. Consensus functional connectivity was introduced here, which was defined as the functional connectivity feature appearing in each cross-validation iteration. Region weight, representing the relative contribution to the identification of schizophrenic patients, was denoted by its occurrence number in the consensus functional connections in this study. The consensus functional connectivity discriminative power was denoted by the mean of its discriminative powers across all iterations of the cross-validation.

Support vector classification and performance evaluation
When the data set of features with high discriminative power was obtained, support vector machines (SVM) with linear kernel function were employed to solve the classification problem (Vapnik, 1995;Bishop, 2006). Due to the classification results were influenced by the number of involved features, classification accuracies with fixed parameter setting of SVM (LIBSVM3.11: Linear kernel, C = 1) using a wide range of feature number were reported. Due to our limited number of samples, we used a leaveone-out cross-validation strategy to estimate the generalization ability of our classifier. The performance of a classifier can be quantified using the generalization rate (GR, i.e., the rate of all the subjects correctly classified), sensitivity (SS, i.e., the rate of the patients correctly classified) and specificity (SC, i.e., the rate of the controls correctly classified) based on the results of cross-validation.

RESULTS
We constructed whole brain functional connectivity using the PCC, MIC, and eMIC based on fMRI data collected from schizophrenic patients and matched healthy controls. Multivariate pattern classification was then introduced to determine the discriminative abilities of the three kinds of functional connectivity. Finally, we analyzed the abnormalities of non-linear and linear functional connections in schizophrenia and determined the spatial distribution of the brain regions related to symptoms of schizophrenia.

CLASSIFICATION
We examined the whole brain resting-sate functional connectivity of the schizophrenic patients and the healthy controls using PCC, MIC, and eMIC, respectively. The 6670 mean functional connections of the controls (Figure 1, the first row) and the patients (Figure 1, the mid row) according to the PCC, MIC, and eMIC. In the third row of Figure 1, the first two panels were presented the relationship between the MIC and PCC, and between the eMIC and PCC, where each red star represented one of the 6670 mean functional connections over the patients, each blue cross represented one of the 6670 mean functional connections over the controls. The last panel in the third row of Figure 1 showed the mean and deviation of the 6670 mean function connections of the patients (red bar) and the controls (blue bar) which were depicted in the first two panels of the third row. Linear SVM was employed to differentiate the patients with schizophrenia from healthy controls using the whole brain functional connectivity extracted by the PCC, MIC, and eMIC.
To evaluate the discriminative power of the PCC, MIC, and eMIC, we conducted multivariate classification. When the numbers of connections used for classification changed, the classification accuracies changed accordingly. The subjects were classified by using the first connections, and 6670 accuracy rates were obtained. With this approach, the eMIC exhibited excellent performance and the highest classification accuracy rate among the three methods. When the number of features for classification is fixed, the obtained accuracy cannot fully reflect the discriminative power of the whole feature space. Here, we present the classification results using a wide range (50, 100, 150, . . ., 1000) of the first selected connections (Figure 2). Clearly, the classification accuracies (GR) of the eMIC were consistently higher than the accuracies of the other two methods across the full feature space.

CHANGES OF FUNCTIONAL CONNECTIVITY
From Figure 2, we can see that when no more than 150 features were involved in the classification, the PCC, MIC, and eMIC gave the best performance. Then, if more features were involved, the classification accuracy of the three  methods decreased and the classification accuracy tend to be stable when the number over 200. Therefore, consensus functional connections of the first 200 features involved in the classification corresponding to the PCC, MIC, and eMIC were evaluated. Additionally, the number of features was identified in accordance with previous study (Dosenbach et al., 2010). Then, 113 consensus functional connections were obtained from each iteration of the leave-one-out cross-validation for the MIC and PCC, and 122 consensus functional connections were identified for the eMIC. The consensus functional connections from the PCC, MIC, and eMIC were projected to a surface rendering of a human brain that was visualized with BrainNetViewer (http://www.nitrc.org/projects/bnv/) (see Figure 3). Further, we added a probabilistic atlas of the Frontiers in Human Neuroscience www.frontiersin.org October 2013 | Volume 7 | Article 702 | 5 FIGURE 3 | Region weights and strength of the connections constructed with PCC, MIC, and eMIC. The connections are displayed in a surface rendering of a human brain. The thicknesses of the consensus connections during leave-one-out cross-validation are scaled by their strength (normalized mean of the rank orders in all iterations during the leave-one-out cross-validation). Connections whose strength were increased in the patients relative to the controls are shown in orange, and connections whose strength were decreased in patients are shown in light green. The ROIs related to the selected consensus connections are also scaled by their weights (sum of the weights of all connections to and from that ROI) and displayed. The ROIs are color-coded according to the functional networks (CON, white; DMN, green; cerebellum, red; visual network, yellow; sensorimotor network, cyan; frontal-parietal network, rose; and other, black). The number in this figure for ROIs is shown in Table A1.
cerebellum to the ICBM152 cerebrum template that was released with the software. Primary consideration was given to changes in the strength of the connections. The consensus functional connections according to PCC and MIC were all decreased in patients with schizophrenia compared to healthy controls. It is generally accepted that the strength of functional connectivity is decreased in patients with schizophrenia compared with healthy controls. The PCCand MIC-based functional connections found in our study are consistent with this common view. In contrast, the strength of the eMIC connections were elevated in patients with schizophrenia. In the Discussion section, we provide a detailed explanation. Additionally, in terms of anatomical location, both the consensus MIC connections and the consensus eMIC connections were largely consistent with the consensus PCC connections. The MIC and the PCC shared 56 common connections, the eMIC and the PCC shared 53 common connections, and the common connections make up approximately half of the total consensus connections. Significant difference of eMIC to PCC is related to the functional connections between the cerebellum and the temporal cortex located in the default model network (DMN) that were found by eMIC but were not detected by PCC. Although there were no PCC functional connections between the cerebellum and the temporal cortex, the cerebellum and the temporal cortex were both connected with the parietal cortex (Figure 4).

DISTRIBUTION OF FUNCTIONAL CONNECTIVITY
In addition, we investigated the distribution of brain regions that are related to the consensus connections produced by the PCC, MIC, and eMIC. Overall, the distributions of ROIs identified by the PCC, MIC, and eMIC were similar, although the consensus connections were not completely identical across the three methods (Figures 3, 4).
To facilitate analysis, these ROIs can be categorized as follows: (i) the cingulo-opercular network (CON), including several regions in the anterior prefrontal cortex, the inferior parietal Frontiers in Human Neuroscience www.frontiersin.org October 2013 | Volume 7 | Article 702 | 6 cortex, the basal ganglion, the dorsal anterior cingulate cortex (dACC), the insular, the thalamus, and the cerebellum (which will be discussed later separately due to its importance), which is a newly defined cognitive network with great significance for schizophrenia (Dosenbach et al., 2007;Tu et al., 2012); (ii) the DMN, including structures of the hippocampus, the posterior cingulate cortex, the medial prefrontal cortex (mPFC), and the bilateral inferior parietal cortex, which is believed to play an important role in the pathogenesis of schizophrenia (Raichle et al., 2001;Fransson, 2005;Whitfield-Gabrieli et al., 2009); (iii) the cerebellum network, which can be seen as part of the CON; (iv) the visual network comprising the primary visual cortex, the extra-striate visual areas and the lingual gyrus, fusiform gyrus, and calcarine gyrus, which is involved in visual processing Damoiseaux et al., 2006;van den Heuvel et al., 2008;van den Heuvel and Hulshoff Pol, 2010a,b); (v) the sensorimotor network, including the primary sensory cortex, the primary motor cortex the supplementary motor cortex (Biswal et al., 1995;van den Heuvel and Hulshoff Pol, 2010a,b); and (vi) the frontal-parietal network, consisting of the superior parietal and the superior frontal cortex, which is involved in attention processing (Dosenbach et al., 2010;De Luca et al., 2006) (Figure 4).
For the three kinds of functional connectivity, the ROIs with the heaviest weight were distributed primarily in the DMN, cerebellum, visual network, and frontal-parietal network. Specifically, the cerebellum was the most important network according to  Figure 3. The green lines represent connections constructed by the PCC, blue lines represent connections constructed by the MIC, and orange lines represent connections constructed by the eMIC. This figure was plotted using the MATLAB toolbox called PlotPie which was developed by our study group and will be released in the near future. eMIC, while the visual cortex was the most weighted network according to PCC and MIC.

DISCUSSION
In this study, we introduced a novel measure called eMIC for estimating the non-linear functional connectivity underlying schizophrenia and applied this estimation of functional connectivity to distinguish schizophrenic patients from healthy controls. Then, we found that strength of the non-linear functional connectivity increased in patients with schizophrenia which was opposed to that observed for the traditional method, which can be attributed to the compensatory mechanisms in the human brain. Furthermore, the non-linear and linear functional connectivity presented similar but not completely the same spatial distribution.

ANALYSIS OF CLASSIFICATION
The results of the classifications produced by the PCC, MIC, and eMIC were obtained using the same procedure and the same classifier parameters. The only factor affecting the classification accuracy is the measure of the functional connectivity. Using support vector classification, we compared the discriminative powers of the three kinds of functional connectivity produced by PCC, MIC, and eMIC. The eMIC produced consistently higher classification accuracies than the other two methods across the full connection space (Figure 2). Furthermore, when all 6670 features were used for classification, we obtained the classification results displayed in Table 2. The classification accuracy of the eMIC remained higher than those of the other methods. In short, the functional connectivity produced by the eMIC had the greatest discriminative ability among all three methods. Now, we pay more attention to the reason for the better classification accuracy of eMIC. The MIC maximizes the association between two time series, whereas the PCC captures only linear function. Then, the eMIC may capture the subtle non-linear neuronal synchronization in human brain, which will improve the performance of eMIC. From another point of view, the eMIC combines the discriminative information of both the PCC and MIC, which reflects a feature-level information fusion that increased classification accuracy. An MI-based study reported that the linear functions accounted for most of the associations between the fMRI time series and that non-linear functions only mined 5% more information in fMRI time series (Hlinka et al., 2011); this could be explain why the eMIC gave limited improvement to PCC (Figure 2). Conservatively, this result confirmed that the non-linear associations have discriminative power which should not be overlooked. We believe that these non-linear connections add new discriminative information to the linear connections which may increase the classification accuracy. Here, we give a detailed analysis of the increased strength of the selected eMIC functional connections. The decreased strength of the PCC functional connections in schizophrenia have been noted in an overwhelming majority of studies (Camchong et al., 2011;Pettersson-Yeo et al., 2011;Repovs et al., 2011;Fornito et al., 2012), which was also confirmed by our previous study (Shen et al., 2010). However, the eMIC-based functional connections demonstrated increased strength in schizophrenic patients compared to healthy controls. As MIC captures various associations but the majority of associations were linear which can be estimated by PCC (Hlinka et al., 2011), decrease of the PCC strength leads to decrease of the MIC strength. The human brain is an organ with great plasticity and adaptability. Thus, the compensatory mechanism of the human brain, that's non-linear correlation will strengthened to reconcile the influence of the corrupt of linear correlation, can result in increases of non-linear functional connectivity. Although the connection selecting procedures of the three methods were implemented separately, the consensus functional connections by the eMIC and the PCC still had a relatively high identity (approximately 50%), which support the notion that the compensatory mechanism of the human brain results in the increase strength of the eMIC in patients with schizophrenia.
If the strength of MIC did not change between the patients and the controls, while the PCC decreased, the increase of the strength of eMIC will be purely caused by the changes of PCC. In fact, the connection strength of MIC decreased in the patients (Figure 1) and gave classification accuracy over 80% (Figure 2). Thus, the increase of eMIC in schizophrenic patients was not simply caused by the decrease of PCC but by the decrease of both PCC and MIC. We believe that the eMIC captures subtle changes of functional connectivity and adds new discriminative information to the classification.

NETWORK ANALYSIS OF NON-LINEAR FUNCTIONAL CONNECTIVITY
As the eMIC provided useful discriminative information to the linear functional connectivity, aberrant non-linear and linear functional connections were both categorized into six networks for explaining the symptom of schizophrenia. Therein, connections between the DMN and the cerebellum, the DMN and the visual network, the DMN and the frontal-parietal network, and the cerebellum and the frontal-parietal network demonstrated the greatest discriminative power for both the linear and non-linear functional connectivity.
The DMN, the CON (cerebellum), visual network, and frontal-parietal network which are related to specific brain functions were all weighted by the linear and non-linear measures. The DMN was frequently reported in previous studies and weighted in our study. The three methods all identified the important nodes in the DMN such as the mPFC, inferior parietal cortex, Para hippocampal gyrus, and middle temporal cortex. This network is generally accepted as an important network associated with schizophrenia (Bluhm et al., 2007;Garrity et al., 2007;Whitfield-Gabrieli et al., 2009). The CON is believed to support the "task model" in opposition to the "default model" (Dosenbach et al., 2006). In our study, the cerebellum_6_R which can be seen as part of the CON is important for all the three functional connectivity, particularly for the new non-linear measure. Additionally, the PCC, MIC, and eMIC all gave the visual network heavy weight. The frontal-parietal network was also identified by the all the three measures, which was reported as important brain regions in previous studies (Honey et al., 2005;Lynall et al., 2010).
The DMN, the CON (including cerebellum), visual network, and frontal-parietal network obtained slightly different region weight according to the PCC, MIC, and eMIC, respectively. Functional connections between the DMN and the cerebellum were identified by the eMIC but not by the PCC. If the cerebellum is viewed as part of the CON, interactions between the DMN and the CON can be established by the eMIC. The wide spread changes in the CON and DMN may interact with each other, and contribute to the functional basis of schizophrenia. In the connections produced by the PCC and MIC, the right inferior occipital cortex in the visual network exhibited the greatest region weights. However, the cerebellum_6_R exhibited the greatest region weight for the eMIC. This result implies that the cerebellum may play a role in non-linear interaction between different brain regions. Connections between the visual network and regions of the frontal cortex are believed to be involved in visual perception which may contribute to the aberrance of visual perception in schizophrenia (Harvey et al., 2011;Calderone et al., 2013). The cerebellum and the DMN are linked to the frontal-parietal network by eMIC and MIC but not PCC. Combined with the DMN, this network is thought to be closely related to attention tasks (Bush and Shin, 2006;Gao and Lin, 2012), which may explain the occurrence of attention impairment and its relationship to the DMN in schizophrenia.
In conclusion, from the functional network perspective, the distribution of ROIs with the greatest weight according to the linear and non-linear connections was similar but not completely the same, and the non-linear connections shed new light on interpretation of the schizophrenic symptoms.

LIMITATIONS
Although the new functional connectivity constructed by using the eMIC method exhibits better performance in the classification and has explored new information about schizophrenia, there are several limitation in this study. First, our sample size was small. Generalizations of the proposed methods need to be tested with large data sets. Second, the number of slices in the fMRI image used in this study was 20, which was relatively fewer and may not sufficiently to capture the details of the abnormalities in the patients. Third, network analysis of the human brain is a trend in the literature which include ROI localization and connectivity estimation. Our study is confined to connectivity analyses that do not include definition of ROIs.

CONCLUSIONS
In this study, we introduced a novel non-linear functional connectivity for schizophrenia study. The classification results show that the non-linear functional connectivity has an equal if not better discriminative ability than existing linear measures in schizophrenia identification. This result suggests that non-linear Frontiers in Human Neuroscience www.frontiersin.org October 2013 | Volume 7 | Article 702 | 8 functional connectivity should be taken into account in research on schizophrenia and other psychiatric disorders. Furthermore, we found that the non-linear functional connectivity which was strengthened in the patients has a similar distribution with its linear counterpart. This new finding indicates the presence of compensatory mechanism between linear and non-linear associations and non-linear functional network abnormalities underlying the symptom of schizophrenia.

STUDY FUNDING
This study was supported by the National Basic Research Program of China (2011CB707802) and the National Natural Science Foundation of China (61003202, 61005084).