A total of 59 right-handed participants, including 20 SIVCIND patients, 20 SIVaD patients, and 19 healthy controls, were recruited. All subjects provided written informed consent in accordance with the Declaration of Helsinki. This study was approved by the local Medical Ethics Committee at the Third Military University (Chongqing, China) on Human Studies.

All of the patients had subcortical WM hyperintensity on T2-weighted imaging with at least two lacunar infarcts, in accordance with the criteria suggested by Galluzzi et al. (17). However, the subcortical WM hyperintensity was not quantified. The Patients with SIVCIND were diagnosed based on the following criteria (18, 19): (1) subjective cognitive complaints reported by the participant or his/her caregiver; (2) insufficient cognitive impairment to meet the DSM-V criteria for dementia; and (3) neuropsychological tests including Hachinski Ischemic Score (HIS) determination (HIS ≥ 7). SIVaD patients were diagnosed on the basis of both DSM-V criteria and the National Institute of Neurological Disorders and Stroke–Association Internationale pour la Recherche et l’Enseignement en Neurosciences (NINDS-AIREN) criteria (20). Exclusion criteria included aphasia and inability to complete psychological testing because of a language disorder, a diagnosed chronic or degenerative disease or condition before stroke affecting the central nervous system, active substance abuse disorders, intracranial hemorrhage and cortical infarct, or epilepsy. However, the exclusion criterion was not including past or present medication for the treatment of cognitive impairment and infarction. Clinical performance was evaluated through the following neuropsychological tasks: Mini-Mental State Examination (MMSE), Clinical Dementia Rating (CDR), Global Deterioration Scale (GDS), Activities of Daily Living Scale (ADL), Montreal Cognitive Assessment (MoCA), and HIS (21–27). All patients received complete sociodemographic and clinical data collection, and two experienced neurologists performed the diagnoses for the two groups.

Nineteen age- and gender-matched healthy volunteers with normal cognitive performance and no known nervous system diseases were recruited as the healthy controls. All healthy controls underwent neurological and neuropsychological evaluation and conventional MRI. None of them had vascular risk factors, cognitive complaints, current psychiatric illness, or a history of psychiatric illness. No healthy controls had brain trauma, brain tumor, psychiatric disorders, systemic disease, or other MRI contraindications such as claustrophobia.

All images were acquired on a 3 T scanner (Magnetom Trio; Siemens Medical Systems, Erlangen, Germany) equipped with eight-channel, phase-array head coils. Foam padding was used to minimize head motion by all subjects. Conventional MRI included transverse T1-weighted images (TR/TE = 200/2.78 ms, flip angle = 70^◦, matrix = 384 × 384, thickness = 4.0 mm, 25 slices, voxel size = 0.7 mm × 0.6 mm × 5 mm) and fluid-attenuated inversion recovery (TR/TE/TI = 9,000/93/2,500 ms, flip angle = 130°, matrix = 256 × 256, thickness = 4.0 mm, 25 slices, voxel size = 0.9 mm × 0.9 mm × 4 mm). Anatomical image datasets were acquired with a magnetization-prepared rapid gradient echo sequence (TR/TE/TI = 1, 900/2.52/900 ms, flip angle = 9◦, matrix = 256 × 256, thickness = 1.0 mm, no gap, 176 slices, voxel size = 1 mm × 1 mm × 1 mm). Resting-state images were acquired in the axial orientation using a single-shot, gradient-recalled echo planar imaging (EPI) sequence using the following parameters: TR = 2,000 ms, TE = 30 ms, flip angle = 90°, matrix size = 64 × 64, FOV = 192 mm × 192 mm, 36 transverse slices, 3-mm slice thickness without a gap, and resolution = 3 mm × 3 mm × 3 mm. During EPI data acquisition, the subjects were instructed to remain awake, relax and keep their eyes closed without thinking of anything in particular. For each subject, we collected 240 volumes in total.

Data preprocessing was performed using Statistical Parametric Mapping software (SPM8, http://www.fil.ion.ucl.ac.uk/spm) and Data Processing Assistant for Resting-State fMRI (DPARSF, http://www.restfmri.net). We discarded the first ten volumes to ensure steady-state longitudinal magnetization, and the remaining 230 volumes were corrected for delay in slice acquisition. The images were later coregistered to the first image for correction of rigid-body head motion. For all subjects, excessive motion was defined as more than ±2 mm and ±2° of the translational and rotational parameters. Subsequently, the individual three-dimensional T1-weighted images were coregistered to the mean functional images, which were generated by SPM in the preprocessing step of alignment, after motion correction using a linear transformation. Additionally, we employed multiple linear regression to remove the signal from white matter and cerebrospinal fluid (28, 29). Then, the functional scans were spatially normalized to the standard Montreal Neurological Institute space and resampled to 3 mm isotropic voxels using an optimum 6-parameter affine transformation and nonlinear deformations. As previously suggested, we did not apply spatial smoothing to avoid introducing artificial local spatial correlations (30, 31). Finally, a band-pass filter (0.01–0.08 Hz) was performed to reduce low-frequency drift and high-frequency physiological noise in each voxel (32).

The nodes and edges between nodes are the two components of a network. In the functional brain network, the nodes represent brain regions. In our study, to determine the nodes of the functional brain network, we used an automated anatomical labeling (AAL) atlas (33) to segment the brain into 90 regions of interest (ROIs) (45 for each hemisphere; details are provided in Table S1 in Supplementary Material). Herein, the functional network was obtained from the defined ROIs and BOLD signal. The nodes represented the 90 ROIs, and the edges represented the Pearson correlation between the mean BOLD time series in two given ROIs. For each participant, a 90 × 90 symmetric weighted network was constructed. Fisher’s r-to-z transformation was applied to the correlation matrix of each subject to improve the normality of the correlation coefficients (34). Finally, as age, gender, and formal education duration differed across groups, the weighted network was corrected for all three by regressing out these factors from the weighted edges of each participant.

The corrected weighted network was thresholded at different levels of sparsity ranging from 10 to 40% using increments of 1%, keeping the highest weighted edges. At each threshold, network metrics were calculated. The functional organization indices of brain networks were computed using the GRETNA toolbox (https://www.nitrc.org/projects/gretna/).

To characterize the weighted brain functional network, several global network parameters were calculated, including the clustering coefficient (C_p), mean path length (L_p), global efficiency (E_glob), and local efficiency (E_loc). The definitions of these metrics are briefly described below.

The clustering coefficient C_p of a network G of N nodes describes the connectedness of direct neighbors around individual nodes (34). The clustering coefficient C_p is expressed as follows (35): Cp(G)=meani∈node2ki(ki−1)∑j,k(wij¯wjk¯wki¯)1/3, where k_i is the degree of node i, (∑j≠i∈Gwij>0), and w is the weighted value of the edge, which is scaled by the mean of all weights to control each participant’s cost at the same level.

The mean path length L_p of G quantifies the ability to propagate information in parallel. The mean path length L_p is expressed as follows: Lp(G)=1N(N−1)∑i≠j∈GLij, where L_ij is the shortest path length between any pair of nodes i and j and is defined as the sum of the edge weight w_ij along this path.

The E_g of G quantifies the global efficiency of parallel information transfer in the network. The E_g of G is expressed as follows: Eg(G)=1N∑i∈GEnodal_glob(i), where Enodal_glob(i)=1N−1∑i≠j∈G1Lij, where L_ij is the shortest path length from node i to j and is defined as the sum of the edge weight w_ij along this path.

The E_loc of G quantifies the fault tolerance of a network, indicating the capacity for information exchange within each subgraph when the index node is eliminated (30). The E_loc of G is expressed as follows: Eloc(G)=1N∑i∈GEnodal_loc(i), where Enodal_loc(i)=Eg(Gi), where G_i denotes the subgraph composed of the nearest neighbors of node i.

The nodal global efficiency E_{nodal_glob} and nodal local efficiency E_{nodal_loc} measure the capacity of information propagation of the given node with all other nodes in network and their direct neighbors, respectively.

Furthermore, we calculated the area under the curve (AUC) for each network metric (C_p, L_p, E_glob, E_loc, E_{nodal_glob}, and E_{nodal_loc}) to provide a summarized scalar for the topological characterization of brain networks independent of single threshold selection (36).

To determine whether a combination of graph-theory measures could predict the subject groups, we performed Adaptive Boosting, a machine learning meta-algorithm used to build classification models. We chose a decision tree as a weak learner, and the number of learning cycles was 100. In our model, the predictors were the four global network metrics (C_p, L_p, E_g, and E_loc), and the response was the subject’s group. Classification accuracy was evaluated via leave-one-out cross-validation to ensure a relatively unbiased estimate of the generalization power of the classifiers to new subjects. The sensitivity and the specificity were calculated from the confusion matrix.

To further identify functional connections showing differences in patient groups, we used the network-based statistics (NBS) approach, which is a validated, nonparametric statistical approach for controlling family-wise error in connectome analyses (37). NBS was applied to compare the patient and normal control groups with the following sequence: (1) one-way ANOVA (post hoc: two-sample one-tailed t-tests) was performed to compare the strength of the edge weight at each individual element of the connectivity matrix, (2) network components of connected edges that survived a p-value of 0.005 uncorrected were retained, and (3) the size of the largest cluster was calculated. To generate an empirical null distribution for evaluating the statistical significance of the cluster sizes, we randomly shuffled the groups (10,000 permutations) and obtained the largest cluster size null distribution by repeating steps 1, 2, and 3. Finally, for any connected component of size M that found in the right grouping of controls and patients, the corrected p-value was determined by calculating the proportion of the 10,000 permutations for which the maximal connected component was larger than M.

Corrected weighted networks of all groups were entered into one-way ANOVA with the NBS approach (p < 0.01 NBS corrected for multiple comparisons). Post hoc two-sample one-tailed t-tests were performed to assess between-group differences in the significant network obtained by ANOVA (p < 0.01 NBS corrected for multiple comparisons).

For visualization purposes, circular graphical representations were used to display significant connections in the statistical analyses using circos software (http://www.circos.ca) (38). Pair-wise connections were displayed with links colored by connection type as follows: left intrahemispheric (blue), right intrahemispheric (green), and interhemispheric (red). ROIs were grouped according to Salvador (i.e., frontal, temporal, parietal, medial temporal, occipital, and subcortical) (39). The ROIs with a high degree of significant connections (k > 1 SD above the mean) were classified as “network hubs” and represented in blue, in accordance with a previous study (10).

All statistical analyses were performed using SPSS Statistics 22.0.0 (http://www-01.ibm.com/software/analytics/spss/).Pearson’s chi-squared test was used to compare sex categorical variables. Three-level one-way ANOVAs were used to compare clinical and sociodemographic data between the NC and two patient subgroups. Two-sample t-test was used to detect differences in clinical scale scores (HIS, GDS, CDR, and ADL) between SIVCIND and SIVaD patients. The statistical significance threshold was set at p < 0.05.

Between-group differences in AUC values of global network metrics (C_p, L_p, E_g, and E_loc) and nodal efficiency (E_{nodal_glob} and E_{nodal_loc}) were investigated with one-way ANOVA after the normality of the data distribution was assessed with a Kolmogorov–Smirnov test. Levene’s test showed nonhomogeneity of variances within groups; thus, Welch ANOVA was performed on the data. Metrics showing a main effect of group differences in the ANOVA model were further evaluated by post hoc tests (Games-Howell test). A significance threshold of p < 0.05 was applied to each test, and the FDR was used to correct for multiple comparisons.

In addition, to explore the relationships between topological properties of the network measures and clinical outcomes in cognitive function, we further performed Spearman’s correlation between topological properties (C_p, L_p, E_g, and E_loc) and cognitive test scores (MMSE and MoCA) in patient groups containing both SIVCIND and SIVaD patients. Z-scores were calculated for each of these cognitive scores. We used a statistical significance level of p < 0.05, and the FDR was used to correct for multiple comparisons.

Subject demographics are listed in Table 1. Significant differences between the three groups (NC, SIVCIND, and SIVaD) were found in MMSE (F = 93.322, p < 0.001) and MoCA scores (F = 115.981, p < 0.001). Games-Howell tests revealed that SIVaD patients and SIVCIND patients showed lower MMSE and MoCA scores than did normal controls, and SIVaD patients showed lower scores than did SIVCIND patients. In addition, two-sample t-test revealed that SIVCIND patients showed lower HIS (T = −3.498, p = 0.001), GDS (T = −5.121, p < 0.001), and CDR (T = −6.658, p < 0.001) scores than did SIVaD patients. There was no significant difference among the three groups in gender (χ² = 0.923, p = 0.632), age (F = 0.663, p = 0.519), and duration of formal education (F = 1.695, p = 0.193).

There was a significant difference between the groups as determined by one-way ANOVA for the clustering coefficient C_p (F = 3.761, p = 0.033), mean path length L_p (F = 9.042, p = 0.001), global efficiency E_g (F = 11.069, p < 0.001), and local efficiency E_loc (F = 10.472, p < 0.001) (Figure 1). A Games-Howell test revealed that E_g and E_loc were significantly lower for SIVCIND and SIVaD than for NC. Additionally, SIVaD showed lower E_g and E_loc values than did SIVCIND. The C_p values were significantly higher for NC than for SIVaD. The L_p values were significantly higher for SIVaD than for SIVCIND and NC.

Compared with NC, SIVCIND patients showed decreased nodal global efficiency in several regions located in the parietal and temporal cortices (Figure 2A). In addition to these regions, SIVaD showed decreased nodal global efficiency in the prefrontal and occipital cortices and more temporal areas (Figure 2A). Moreover, compared with SIVCIND patients, SIVaD patients exhibited lower nodal global efficiency in the right precentral gyrus, cuneus, superior temporal gyrus, and inferior temporal gyrus (Figure 2A).

Additionally, compared with NC, SIVCIND patients showed decreased nodal local efficiency predominantly in the right putamen, superior parietal gyrus, supramarginal gyrus, angular gyrus, and bilateral amygdala (Figure 2B). SIVaD had widely decreased nodal local efficiency in the whole brain (Figure 2B). Compared with SIVaD patients, SIVCIND patients had higher nodal local efficiency in the right precentral and postcentral gyrus and some default-mode regions (Figure 2B). Compared with NC, patients did not show a significant increase in nodal global efficiency and nodal local efficiency.

ANOVA revealed significant functional differences between groups (NBS corrected, p < 0.01) in almost the whole brain (Figure 3; Tables 2 and 3). Post hoc analysis revealed between-group differences in the architectural features of connections and hub nodes; these differences are illustrated by connectograms (Figure 4; Tables 2 and 3). Compared with NC, SIVCIND showed 93 significantly decreased connections, and SIVaD exhibited 141 significantly decreased connections. Notably, the temporal and medial temporal cortices played more important roles in significantly decreased connections in SIVaD (39%) than in SIVCIND (26.8%). Additionally, compared with SIVCIND, SIVaD exhibited 34 significantly decreased connections, which were mostly located in the frontal/temporal and temporal/temporal cortices.

Moreover, compared with NC, SIVCIND had lower hub connections in the precuneus and supramarginal gyri in the bilateral hemisphere; in the middle cingulate gyrus, rolandic operculum, and superior temporal gyrus in the right hemisphere; and in the inferior parietal gyrus and orbital middle frontal gyrus in the left hemisphere. SIVaD had lower hub connections in the superior temporal gyrus and supramarginal gyri in the bilateral hemisphere and in the paracentral lobule, middle temporal gyrus, superior frontal gyrus (medial orbital), and middle occipital gyrus in the right hemisphere. In addition, compared with SIVaD, SIVCIND had higher hub connections in the superior frontal gyrus (medial orbital and orbital areas) and some temporal areas.

Correlation analysis showed that MMSE and MoCA scores were positively correlated with C_p, E_loc, and E_g values and negatively correlated with L_p values (Figures 5 and 6).

The leave-one-out cross-validation classification procedure was applied using the four network metrics, and the performance of the procedure is summarized in Table 4. A classification accuracy of 76% was obtained for all subjects. However, the accuracy varied between groups. The procedure was good for distinguishing SIVaD patients from normal controls with a sensitivity of 95% and a specificity of 90%. More difficulties were observed when distinguishing SIVCIND patients from normal controls with a sensitivity of 89% and a specificity of 65%.

As shown by the progressive decrease in the values of the clustering coefficient, local efficiency, and global efficiency as cognitive impairment increased, there was a reduction in the functional integration and functional segregation of the brain with cognitive impairment. Functional segregation in the brain represents specialized processing within densely interconnected groups of brain regions, and functional integration in the brain represents the ability to rapidly combine specialized information from distributed brain regions (40). This information suggests that progressive cognitive impairments in SIVD may be related to the progressive decrease in the density of connections involved in specialized networks engaged in cognition and the progressive decrease in long-distance functional connections. Specifically, MMSE and MoCA scores were positively correlated with functional integration and segregation of the brain. Hence, decreased performance on the cognitive tests among patients with cognitive impairment caused by SIVD was associated with less extended local clustering, less efficient communication among neighboring areas, and a reduced ability to combine specialized information from distributed brain regions. This finding is in line with the results of previous studies investigating functional organization. For instance, Yu et al. (16) found that compared with normal controls, subcortical vascular cognitive impairment patients exhibited a decreased clustering coefficient and global efficiency in parietal and temporal areas; these changes were correlated with cognitive test performance. By selecting the thalamus as the seed, Zhou et al. (9) reported that MMSE scores correlated with decreased functional connectivity between the left thalamus and left olfactory bulb in subcortical vascular cognitive impairment patients.

In addition to global metrics, we also studied local and global efficiency at a nodal level to show the extent of information transmission capacity of nodes with their neighbors and all other nodes in networks. As shown in Figure 2, the brain regions with altered nodal efficiency progressively extended over the whole brain as cognitive impairment worsened. Compared with NC, SIVCIND and SIVaD did not exhibit significantly increased nodal efficiency, which suggested reduced information transmission between a node and its neighbors and all other nodes in the network in patients with cognitive impairment caused by SIVD. Specifically, we found that SIVCIND showed decreased nodal efficiency in several brain regions involved in the parietal/temporal cortex. When cognitive impairment increased (SIVaD), more brain regions with decreased nodal efficiency were observed in the parietal/temporal lobe. Structural abnormality evidence has shown that the cortex thins in the parietal and temporal cortices as cognitive impairment worsens (4). Thus, this finding is compatible with previous studies, and the gradual progression of decreased nodal efficiency in parietal/temporal areas suggested that their lower capacity for information exchange with other nodes in the network presumably occurred in response to non-cognitive consequences such as apathy and slower motor function in patients with cognitive impairment caused by SIVD (3, 41, 42).

Moreover, SIVCIND had lower nodal global efficiency in the lateral DMN, such as the precuneus and left posterior cingulate gyrus, in line with previous studies (9, 43). When cognitive impairment increased (SIVaD), the regions of significantly decreased nodal global efficiency extended to the prefrontal cortex, including the medial superior frontal gyrus and orbital and dorsolateral frontal cortices, which are the primary components of the prefrontal/subcortical circuit, and the interruption of this circuit was associated with executive dysfunction (1, 3), which is a specific cognitive profile in cognitive deficits affected by vascular lesions (44). Additionally, a similar pattern of progression was observed in the prefrontal cortex with regard to nodal local efficiency as cognitive impairment increased. Thus, the nodal efficiency dysfunction in the prefrontal cortex could be used to determine whether patients are at risk of major cognitive impairment or dementia. Compared with SIVCIND, SIVaD had lower nodal global and local efficiency in the temporal cortex and several occipital areas, such as the cuneus and lingual gyrus; both structural and functional abnormalities in these areas have been associated with depression (14, 45). Depression is particularly prominent in those with vascular dementia (2). Our findings are consistent with these previous studies and suggested decreased information exchange in occipital areas, which presumably occurred in response to depression frequently caused by SIVaD. Unfortunately, the present study included no neuropsychiatric test to estimate the state of neuropsychiatric symptoms such as depression.

The NBS analyses showed that alterations in brain network functional connectivity progressively increased as cognitive impairment worsened (Figure 4). A structural study found that vertices with cortical thinning largely overlapped between SIVaD and SIVCIND, but their extent and severity were greater in SIVaD than in SIVCIND (5, 46). The progressive increase in reduced connectivity with worse cognitive impairment was in accordance with the extent and severity of cortical thinning.

Compared with normal controls, SIVCIND exhibited weakened connections mainly located in the frontal cortex and parietal areas. An examination of the connectogram (Figure 4, top-left) shows that these alterations mainly concerned interhemispheric connections to the frontal and parietal cortices, which were also found in an electroencephalography study (47). Changes in the frontoparietal attentional network (48) agree with the attention deficits usually exhibited by patients with cognitive impairment caused by SIVD. Notably, significantly decreased functional connectivity in the orbitofrontal cortex in SIVCIND is known to be involved in executive function (49, 50). When cognitive deficits worsened (SIVaD), the functional network was widely altered (Figure 4, top-right). Specifically, SIVaD exhibited more decreased connections in the orbitofrontal cortex and dorsolateral prefrontal cortex than did normal controls. The orbitofrontal cortex and dorsolateral prefrontal cortex are the primary components of the prefrontal/subcortical circuit, and the interruption of this circuit was associated with executive dysfunction (1, 3). The progressively increased changes in functional connectivity in the prefrontal cortex were in accordance with the progression in nodal properties and further suggested that the dysfunction of the prefrontal cortex could be used to determine the severity of cognitive decline in SIVD. Additionally, more significantly decreased connections were observed in the temporal and medial temporal cortices for SIVaD than for SIVCIND. Previous studies suggested that the temporal cortex is important in regulating complex social interactions (51) and serves as an important component of the social perception system (52). Therefore, decreased functional connections associated with temporal areas may explain the impaired social behaviors and loss of empathy that are frequently observed in SIVaD. Moreover, the medial temporal cortex is thought to be involved in encoding declarative memory (53), and decreased connections in this area can lead to forgetfulness, a clinical manifestation of SIVaD (3).

Compared with SIVCIND, SIVaD showed weakened connections that mainly concerned interhemispheric connections from the temporal cortex to the temporal, occipital and orbitofrontal cortices (Figure 4, bottom). Additionally, almost all the network hubs were located in temporal areas (Table 3). Particularly, Seo et al. (5) found that SIVaD showed cortical thinning in all these areas involved in SIVCIND in addition to the temporal cortices. Our findings are consistent with previous studies and further suggested that temporal areas may be associated with cognitive decline from SIVCIND to SIVaD.

Several limitations in our study need to be further addressed. First, the functional brain network was constructed at a regional level based on a previously published atlas (AAL) by parcellating the entire brain into 90 regions. Brain networks derived using different parcellation schemes exhibit distinct topological organization (54). In future, different brain parcellation schemes should be used to determine which is most appropriate for characterization of topological architecture in cognitive deficits affected by SIVD. Second, in our study, no neuropsychiatric scale was used to verify the relationship between the abnormal topological properties and neuropsychiatric symptoms observed in patients. However, neuropsychiatric symptoms such as depression are particularly prominent in SIVaD patients. In the future, we will collect detailed data using a neuropsychiatric scale to evaluate the neuropsychiatric impairment associated with subcortical vascular lesions. Third, the sample size in the current study was not large. Although our sample has the important advantage of comprising different levels of cognitive deficits in patients affected by SIVD, our study may have lacked sufficient statistical power to detect some subtle effects among the three groups. Fourth, the number and location of the lesions is not uniform across the patients in our study. In the further studies, more patients will be recruited to ensure the lesions consistently across all patients as far as possible. Finally, similar to most rs-fMRI studies, the interference of certain potential confounding factors, such as respiratory and cardiac cycle artifacts, could not be eliminated by a slow sampling rate, despite using a band-pass filter of 0.01–0.08 Hz.