Resting-State fMRI in Chronic Patients with Disorders of Consciousness: The Role of Lower-Order Networks for Clinical Assessment

Resting-state fMRI (rs-fMRI) is a widely used technique to investigate the residual brain functions of patients with Disorders of Consciousness (DoC). Nonetheless, it is unclear how the networks that are more associated with primary functions, such as the sensory–motor, medial/lateral visual and auditory networks, contribute to clinical assessment. In this study, we examined the rs-fMRI lower-order networks alongside their structural MRI data to clarify the corresponding association with clinical assessment. We studied 109 chronic patients with DoC and emerged from DoC with structural MRI and rs-fMRI: 65 in vegetative state/unresponsive wakefulness state (VS/UWS), 34 in minimally conscious state (MCS) and 10 with severe disability. rs-fMRI data were analyzed with independent component analyses and seed-based analyses, in relation to structural MRI and clinical data. The results showed that VS/UWS had fewer networks than MCS patients and the rs-fMRI activity in each network was decreased. Visual networks were correlated to the clinical status, and in cases where no clinical response occurred, rs-fMRI indicated distinctive networks conveying information in a similar way to other techniques. The information provided by single networks was limited, whereas the four networks together yielded better classification results, particularly when the model included rs-fMRI and structural MRI data (AUC = 0.80). Both quantitative and qualitative rs-fMRI analyses yielded converging results; vascular etiology might confound the results, and disease duration generally reduced the number of networks observed. The lower-order rs-fMRI networks could be used clinically to support and corroborate visual function assessments in DoC.


Introduction
Disorders of Consciousness (DoC) include a spectrum of conditions ranging from the coma, vegetative state/unresponsive wakefulness state (VS/UWS) to a minimally conscious state (MCS), which results from alterations in arousal and/or awareness caused by severe brain damage.
While some patients recover consciousness, others survive in a chronic VS/UWS or MCS condition, with MCS having a better prognosis than VS/UWS [1,2] and the clinical outcome is better for traumatic brain injury than vascular and anoxic etiologies [3][4][5].
The clinical examination based on standardized behavioral scales remains the gold standard for the diagnosis of DoC [6], although two recent guidelines, the American Academy of Neurology [7] and the European Academy of Neurology (EAN) [8], emphasize the role of neuroradiological findings to improve diagnostic accuracy and classification between VS/UWS and MCS. In fact, neuroimaging techniques represent a crucial part of the assessment, providing information about the structural and functional integrity of the patient's brain [9][10][11]. In particular, resting-state functional MRI (rs-fMRI) is applied to the study of DoC to detect residual brain activity. Emerging evidence leads us to consider the usefulness of rs-fMRI in addition to MRI examination as part of a multimodal assessment, as also recommended in the 2020 EAN guideline [8].
The rs-fMRI technique measures the spontaneous neuronal activity that is generated by the brain at rest, reflected by low-frequency fluctuations in the BOLD signal [12]. Regions that show strong synchrony over time are functionally connected and form separate rs-fMRI networks [13]. Several neural networks have been identified and broadly classified into two groups [14]: higher-order networks such as the default mode network (DMN) associated with introspective processes and self-awareness [15], and lower-order networks such as the sensorimotor (SM), visual, and auditory networks (AUD), which exhibit patterns of activity that resemble activities observed during tasks, and are associated with more specialized functions. The SM network, for example, involves regions engaged in sensorimotor functions [16], and its activity at rest displays a degree of hemispheric lateralization that correlates with that one observed during a finger-tapping task [17]. The visual network can be divided into lateral (LVIS) and medial (MVIS) networks, with the first including mesial visual areas, and the latter including lateral occipito-temporal regions. The AUD network consists primarily of Heschl's gyrus and the superior temporal gyrus. In addition, functional connectivity observed at rest in healthy subjects is influenced by experience and can be modulated by prior history of co-activation during active behavior [18][19][20][21]; therefore, functional connectivity in patients with DoC can provide information on the residual functioning of sensory modalities such as the visual, auditory, and somatosensory systems.
In fact, rs-fMRI has been applied to DoC to assess residual neurofunctional activity. The DMN activity was shown to be absent in brain death and generally stronger in MCS than in VS/UWS [22][23][24]. Beyond the well-studied DMN, other rs-fMRI networks might be able to define consciousness, although the number of studies is limited. Among the networks associated with primary functions (SM, LVIS, MVIS, AUD), one study examined 51 cases (26 MCS, 19 VS/UWS and 5 in coma) generally in the subacute phase and validated the classification on another 22 patients with DoC [25]. Each network could discriminate MCS from VS/UWS with high accuracy (>80%), with the AUD having the best performance. Another study involving 16 patients with DoC in the subacute phase also revealed impairment in rs-fMRI connectivity in the SM and AUD networks, but not in the VIS network relative to healthy controls [14]. Lastly, functional connectivity has been studied both within and between networks: rs-fMRI measurements were associated with 1-year clinical outcome in comatose patients [26] and could be used, along with clinical data, to predict prognosis in the chronic phase [27]. Although rs-fMRI can be useful as a supplement to structural MRI as part of a multimodal assessment in DoC [8,28], a wellpreserved map of the rs-fMRI network does not provide evidence of consciousness and the difference between VS/UWS and MCS can be undetected [23,29]. In this context, structural MRI combined with rating scales of anatomical damage has been less investigated in the study of DoC, although recent evidence shows that this type of MRI-based assessment can be useful in discriminating VS/UWS from MCS [24,30].
In regard to the differences in the etiology, a recent multicentric study has shown that post-anoxic patients have the highest clinical complexity, as they have a lower level of consciousness, higher functional disabilities and a greater need for medical devices. This is followed by vascular patients who demonstrate a greater premorbid clinical comorbidity and then by traumatic patients who demonstrate a smaller clinical complexity [31].
However, there is a lack of characterization of these lower-order rs-fMRI networks in terms of their occurrence, etiology, and relationship to behavioral responses. Further, the diagnostic accuracy of networks between VS/UWS and MCS needs to be verified in a larger sample.
In this study, we investigated the integrity of the 4 lower-order rs-fMRI networks (SM, LVIS, MVIS, AUD) associated with primary functions in 109 chronic patients with DoC and emerged from DoC with different etiologies and on average a disease duration > 2 years.
The aim is to better characterize patients in terms of low-order rs-fMRI networks and explore patients' brain residual functional activity in relation to the clinical data, i.e., etiology, disease duration and CRS-R subscores. In addition, the 4 rs-fMRI networks were evaluated alongside their structural MRI, as occurs in clinical practice, and analyzed (i) individually and (ii) by combining the information of the four networks together. We used seed-based analysis (SBA) with map selection among different seeds [32] and independent component analysis (ICA) with the integration of three pieces of information (spatial maps, time series and power spectra) [33,34], developing a qualitative rs-fMRI rating and a quantitative rs-fMRI map intensity index.

Participants
A group of 122 adult patients with DoC and emerged from DoC in chronic phase, admitted to a 1-week program of multidisciplinary assessment at the Coma Research Center (CRC) of the Fondazione IRCCS Istituto Neurologico "Carlo Besta", Milan, Italy, was enrolled. Thirteen patients were excluded for excessive head movement during the MRI session (3) or excessive artefactual MRI data (10). Thus, according to the Aspen criteria, the study included 65 patients in VS/UWS, 34 in MCS and 10 with severe disability (SD) emerged from MCS, as assessed with the CRS-R. Etiology included 33 traumatic brain injury, 39 vascular brain injury and 37 anoxic brain injury. Median disease duration was 27 months (range 2-252, >12 months for 82 cases), and the median age was 50 years (range 19-83). Demographic and clinical characteristics are reported in Table 1. Patients were clinically assessed with the Italian version of the Coma Recovery Scale-Revised CRS-R) [35,36] and with the CRS-R-Modified score [37]; each patient was independently assessed 4 times by experienced raters and the best response was used to establish the final score. Patients also underwent a multimodal assessment comprising evoked potentials (EPs) to assess the presence/absence of the visual, auditory and sensory pathways (as reported in [38]) and FDG-PET, to have a measure of metabolism (as reported in [24,39,40]).
A total of 34 healthy participants (median age 39 years, range 17-66) with no history of neurological deficits were recruited as controls. The study was approved by the ethics committee of the "Carlo Besta" Institute. Written informed consent was obtained from the legally authorized representative of the patients and healthy participants.

MRI Data Acquisition
Scanning was performed on a 3T scanner equipped with a 32-channel head coil (Achieva TX; Philips Healthcare, Best, Amsterdam, The Netherlands). For rs-fMRI data, gradient echo-planar images (EPI) were acquired (repetition time = 2.8 s, echo time = 30 milliseconds, flip angle = 70 • , 2.5 mm isotropic voxel size, 90 × 95 matrix size, 50 slices with 10% gap, ascending order, 200 volumes); the sequence duration was~9.5 min. Sagittal 3D T1-weighted, 2D T2, T2* and FLAIR-weighted images were also acquired. When patient posture allowed, the head was restrained using foam pillows, and a knee wedge was positioned to minimize spine movement and discomfort. Sedation was never performed.

rs-fMRI Data Preprocessing
Data were preprocessed using SPM12, FSL, and in-house code running under Matlab (MathWorks, Natick, MA, USA) [41]. Preprocessing consisted of rigid-body realignment, slice time correction and identification of outlier scans for the scrubbing procedure using a framewise displacement greater than 2 mm (see [25] for similar choices) in FSL volume outlier function. The 3D T1-weighted image was co-registered to mean EPI image, and tissue probability maps of white matter (WM), cerebrospinal fluid (CSF) and grey matter (GM) were extracted.
Multiple linear regressions were used to remove additional confounds, including the 6 head movement parameters, mean WM and CSF signals, as well as linear trend obtained with the 4rd-order polynomial and outlier volumes. A 0.01-0.1 Hz band-pass filtering was applied. The 3D T1 and EPI images were normalized to MNI template. EPI images were spatially smoothed using 6 mm fullwidth at half-maximum Gaussian filter.

Independent Component Analysis (ICA)
ICA was performed using MELODIC (FSL tool https://fsl.fmrib.ox.ac.uk/fsl/fslwiki/ MELODIC, accessed on 31 December 2021) with a fixed number of 30 components [42], which were thresholded at z ≥ 3. Group ICA maps were computed over healthy subjects, and binary templates of the 4 networks were created. Cross-correlation between template and patient single component was computed for each network with fslcc to obtain a similarity index and have a preliminary selection of the more likely neural networks. Then, the component corresponding to each network was identified upon agreement of two experienced observers who considered the spatial map, time series shape, and power spectral density [34]. A component was deemed a candidate network if it exhibited focal activity in the network nodes and specificity to the rest of the brain [43], if the power spectrum presented low frequency only, and if the time series was smooth. Melodic outputs were resampled with a 2 mm voxel.

Seed-Based Analysis (SBA)
For each patient, SBA was performed with DPARSF-A [44] using a 6 mm radius sphere for each node of the 4 networks. A total of 10 seeds was defined, where SM had 4 seeds, 2 located in the left and right pre-and post-central gyri and 2 in left and right supplementary motor area (SMA); AUD had 2 seeds located in the left and right Heschl's gyrus; LVIS had 2 seeds located in left and right inferior occipital and fusiform gyrus; and, finally, MVIS had 2 seeds located in the left and right lingual gyrus and calcarine cortex. The mean time-course of each seed was entered into a whole-brain correlation analysis, and corresponding statistical maps were thresholded at z > 0. 55. In order to select the best seed-based maps, two methods were used to select the location of the seed: the first one was based on the peak coordinates taken from the literature [25,27,45], while the second one was based on subject-specific coordinates placed on the most preserved GM area visible on T1 scan. For each node, the two maps were compared and the less noisy and more specific one was chosen [32,46]. To assess the contribution of each hemisphere and avoid overestimation of connectivity maps, left and right seeds were analyzed separately [47].

rs-fMRI Rating (Qualitative Analysis)
For each patient, a rating was performed by two expert investigators blind to patient data. Each node of ICA and SBA maps was assessed according to the following scale: 0 (no map or undefined map), 1 (map with low specificity and ambiguous cluster), 2 (neuronal map with specific cluster on a few slices) and 3 (neuronal map with highly specific cluster on many slices). For SBA maps the threshold was changed to reduce noise when necessary. In case of large disagreement between the two raters, the scores were reconsidered. The scores obtained by two raters were averaged together. Inter-rater agreement was p = 0.79, indicating good reliability [48].
The presence of a network was defined as a score ≥ 2 on at least one node of ICA or SBA for MVIS, LVIS and AUD networks. For the detection of the SM network, a score ≥ 2 on at least two nodes of ICA or SBA maps was necessary. For each diagnostic group, the number of networks observed was calculated.

rs-fMRI Map Intensity (Quantitative Analysis)
In order to have a quantitative analysis in addition to the rs-fMRI rating, the mean intensity of ICA and SBA maps was extracted for each node, added together and included: (i) as two predictors (left and right hemisphere) in the multivariate logistic regression analysis; (ii) as single variable per network with the Mann-Whitney U test.

Voxelwise Group Map Analyses
Group-level analyses were performed using Statistical NonParametric Mapping (http: //www.nisox.org/Software/SnPM13/, accessed on 31 December 2021) for each network. This was chosen as data were not normally distributed due to the presence of zero-filled maps. The analysis was conducted on SBA maps only because ICA provided fewer recognizable networks. The left and right maps were combined considering for each voxel the maximum value [24]. A correlation analysis was performed between each network and its corresponding CRS-R subscale. Maps with rs-fMRI ratings = 0 or 0.5 were not included to avoid spurious correlations.

MRI Rating
Two expert neuroradiologists blind to all patient data rated independently the severity of gross anatomical and signal abnormality of the regions corresponding to the network nodes, according to a 4-score scale ranging from 0 (severely damaged, i.e., parenchyma obliterated and/or intense, pervasive hyperintensity) to 4 (normal-appearing) (see [24]). Ratings were reconsidered in cases of large disagreement; the scores of each node were averaged between raters and then summed for each network. Inter-rater agreement was p = 0.77, indicating good reliability [48].

Statistical Analyses
Statistical analyses were performed using R software version 4.0.3 (https://cran.rproject.org/, accessed on 31 December 2021) [49]. The Mann-Whitney U test was used to test VS/UWS vs. MCS differences. Spearman ρ was used to study the correlation of imaging data with clinical data, i.e., with CRS-R total score and the CRS-R subscale scores [36], with the CRS-R-Modified score [37] and disease duration. As regards the correlation with CRS-R subscales, the auditory, visual and motor function subscales were selected as they may most correspond to the low-order rs-fMRI networks.
A multivariate logistic regression model was performed to assess the ability to classify VS/UWS vs. MCS patients for single networks and the 4 networks together. The following 5 imaging models were created, with ICA and SBA scores considered separately for rs-fMRI and using 2 scores per node: (1) rs-fMRI rating (2 scores per node); (2) rs-fMRI map intensity (2 scores per node); (3) MRI rating (1 score per node); (4) rs-fMRI rating and MRI rating (3 scores per node); and (5) rs-fMRI map intensity and MRI rating (3 scores per node). Moreover, clinical variables (etiology, disease duration and age) were added in each of the previous models to create imaging + clinical variable models. A model with only clinical variables was used as a baseline to compare all previous models.
The least absolute shrinkage and selection operator (LASSO) method was used to reduce the dimension of explanatory variables of the logistic regression and make the model easier to interpret. Leave-one-out cross-validation (LOOCV) and 10-fold crossvalidation (CV) were used to internally validate the selected model; the results of the 10-fold CV are reported in the Supplementary Material. Accuracy and balanced accuracy obtained after LOOCV were extracted to measure the diagnostic discrimination, together with the estimate of the area under a receiver operating characteristics curve (AUC). The sensitivity and specificity of each AUC were calculated with corresponding CI using the bootstrap method with 2000 stratified replicates. Results are reported in the Supplementary Tables. Following the STARD guidelines [50], false positive and false negative patients were also identified, and diagnostic accuracy was considered as very good for 0.8 < AUC < 0.9, good for 0.7 < AUC < 0.8, sufficient for 0.6 < AUC < 0.7 and bad for 0.5 < AUC <0.6 [51]. The McNemar test was used to compare different logistic regression models for the 4-network assessment.

Differences between Diagnostic Groups
Considering the number of networks observed in each diagnostic group, the most frequent values were 0 for VS/UWS, 2 for MCS and 4 for SD ( Figure 1). Each diagnostic group had a range of variation between 0 and 4 except for SD which never had 0. map intensity (2 scores per node); (3) MRI rating (1 score per node MRI rating (3 scores per node); and (5) rs-fMRI map intensity and node). Moreover, clinical variables (etiology, disease duration a each of the previous models to create imaging + clinical variable only clinical variables was used as a baseline to compare all previ The least absolute shrinkage and selection operator (LASSO reduce the dimension of explanatory variables of the logistic re model easier to interpret. Leave-one-out cross-validation (LOO validation (CV) were used to internally validate the selected mod fold CV are reported in the Supplementary Material. Accuracy obtained after LOOCV were extracted to measure the diagnostic d with the estimate of the area under a receiver operating character sensitivity and specificity of each AUC were calculated with corr bootstrap method with 2000 stratified replicates. Results Supplementary Tables. Following the STARD guidelines [50], f negative patients were also identified, and diagnostic accuracy w good for 0.8 < AUC < 0.9, good for 0.7 < AUC < 0.8, sufficient for 0 for 0.5 < AUC <0.6 [51]. The McNemar test was used to com regression models for the 4-network assessment.

Differences between Diagnostic Groups
Considering the number of networks observed in each diag frequent values were 0 for VS/UWS, 2 for MCS and 4 for SD (Fig group had a range of variation between 0 and 4 except for SD wh Figure 1. Percentage of patients having 0, 1, 2, 3 or 4 networks for each d detection was defined as a score ≥2 on at least one node of ICA or SBA m AUD networks, and defined as a score ≥ 2 on at least two nodes of ICA Network. Figure 1. Percentage of patients having 0, 1, 2, 3 or 4 networks for each diagnostic group. Network detection was defined as a score ≥ 2 on at least one node of ICA or SBA maps for LVIS, MVIS, and AUD networks, and defined as a score ≥ 2 on at least two nodes of ICA or SBA maps for the SM Network. For each network, the Mann-Whitney test showed a significant difference between VS/UWS and MCS (p < 0.05), which was detected more frequently with the rs-fMRI rating than rs-fMRI map intensity ( Figure 2). Correlation with CRS-R total score was significant only with the rs-fMRI rating for the MVIS network (ρ = 0.26 p < 0.01). Considering only VS/UWS and MCS, the correlation was not significant with CRS-R total score, but was significant with CRS-R -Modified score (ρ = 0.22 p < 0.05). Correlation with disease duration was never significant. For each network, the Mann-Whitney test showed a significant difference between VS/UWS and MCS (p < 0.05), which was detected more frequently with the rs-fMRI rating than rs-fMRI map intensity ( Figure 2). Correlation with CRS-R total score was significant only with the rs-fMRI rating for the MVIS network (ρ = 0.26 p < 0.01). Considering only VS/UWS and MCS, the correlation was not significant with CRS-R total score, but was significant with CRS-R -Modified score (ρ = 0.22 p < 0.05). Correlation with disease duration was never significant. For each network, on the left, the spatial map generated from the group-level ICA on control subjects; on the right, group differences measured with rs-fMRI rating and rs-fMRI map intensity between VS/UWS, MCS and SD patients. Boxplot with medians and interquartile range are reported, with Mann-Whitney Z scores. For each network the correlation of rs-fMRI rating with CRS-R total score is reported. * p < 0.05; ** p < 0.01.

Diagnostic Accuracy of Single Networks
Considering the networks individually, logistic regression after applying LASSO showed averaged accuracy around 0.68 with AUC = 0.55 for the rs-fMRI rating and rs-fMRI map intensity models, while balanced accuracy dropped to 0.57 (Table 2,  Supplementary Tables S1 and S2). Generally, there were no differences between the two rs-fMRI models. Combining rs-fMRI with MRI data, accuracy scores improved (on average balanced accuracy = 0.58 and AUC = 0.67) albeit not significantly, and adding clinical variables the results improved slightly further (on average balanced accuracy = 0.64 and AUC = 0.74). The model with only clinical variables (etiology, disease duration and age) has accuracy = 0.72, balance accuracy = 0.64 and AUC = 0.71. Henceforth balanced accuracy is used instead of accuracy because it is considered more appropriate for the presence of unbalanced groups. Finally, there were no hemispherical differences in the number of variables selected by LASSO (see Supplementary Table S3 for the description of the variables). For each network, on the left, the spatial map generated from the group-level ICA on control subjects; on the right, group differences measured with rs-fMRI rating and rs-fMRI map intensity between VS/UWS, MCS and SD patients. Boxplot with medians and interquartile range are reported, with Mann-Whitney Z scores. For each network the correlation of rs-fMRI rating with CRS-R total score is reported. * p < 0.05; ** p < 0.01.

Diagnostic Accuracy of Single Networks
Considering the networks individually, logistic regression after applying LASSO showed averaged accuracy around 0.68 with AUC = 0.55 for the rs-fMRI rating and rs-fMRI map intensity models, while balanced accuracy dropped to 0.57 (Table 2, Supplementary  Tables S1 and S2). Generally, there were no differences between the two rs-fMRI models. Combining rs-fMRI with MRI data, accuracy scores improved (on average balanced accuracy = 0.58 and AUC = 0.67) albeit not significantly, and adding clinical variables the results improved slightly further (on average balanced accuracy = 0.64 and AUC = 0.74). The model with only clinical variables (etiology, disease duration and age) has accuracy = 0.72, balance accuracy = 0.64 and AUC = 0.71. Henceforth balanced accuracy is used instead of accuracy because it is considered more appropriate for the presence of unbalanced groups. Finally, there were no hemispherical differences in the number of variables selected by LASSO (see Supplementary Table S3 for the description of the variables).

Correspondence between Single Networks and CRS-R Subscales
As reported in Table 3, the scores of the CRS-R subscales were centered on score 2 for the motor function and between scores 1 and 2 for the auditory function: this polarization reduced the possibility of observing significant correlations with rs-fMRI. In contrast, the scores of the visual function were more distributed: the visual function was positively correlated with the rs-fMRI rating of the LVIS (ρ = 0.27, p < 0.01) and MVIS (ρ = 0.38, p < 0.001) networks, indicating that the higher the visual score, the more recognizable is the rs-fMRI visual map. A voxelwise group analysis performed between the visual networks (LVIS and MVIS) measured with SBA (which included more maps than ICA) and the CRS-R visual function subscale showed significant results. In particular, for the MVIS network on a total of 88 patients (44 VS/UWS, 34 MCS and 10 SD) results displayed the involvement of the lingual gyrus and calcarine cortex, bilaterally at two thresholds, while for the LVIS network on a total of 84 patients (55 VS/UWS, 23 MCS and 6 SD) results displayed the involvement of the inferior occipital areas more on the left side ( Figure 3). Table 3. Frequency distribution of CRS-R subscale function scores for the different diagnostic categories. Analyzing the correspondence between the CRS-R subscale and network recognizability, we expected an absence of networks when the CRS-R score was zero. However, some networks were detected even when patients had 0 at the subscale: three networks for AUD, six for LVIS and four for MVIS, where the four patients with MVIS had also LVIS ( Table 4). As a result, a total of nine patients had a score ≥ 2 for rs-fMRI rating, indicating high network preservation even if the CRS-R subscale score was 0. Since our patients underwent also visual, auditory and somatosensory potentials-and for some cases also FDG-PET-we compared the results of these techniques in terms of presence/absence of EP and FDG-PET hypometabolism with rs-fMRI data: for seven out of nine cases, the EPs showed the integrity of the pathways and FDG-PET showed active areas in agreement with rs-fMRI results. Thus, concordance was reached for seven out of nine cases, indicating that rs-fMRI data were confirmed by EPs or FDG-PET or both. lingual gyrus and calcarine cortex, bilaterally at two thresh network on a total of 84 patients (55 VS/UWS, 23 MCS and 6 involvement of the inferior occipital areas more on the left side   Table 4. Number of patients with zero score at the CRS-R subscales, for whom rs-fMRI networks was deemed present or absent.

Diagnostic Accuracy of the 4 Networks
Considering the 4 networks, logistic regression after applying LASSO showed averaged balanced accuracy of 0.65 and 0.67 for the two rs-fMRI models ( Table 5, Supplementary  Tables S4 and S5).
Combining rs-fMRI with MRI data, accuracy and balanced accuracy did not change.
Adding clinical variables, accuracy slightly improved to 0.82 for rs-fMRI rating + MRI rating (with balanced accuracy = 0.78) and remained at 0.74 for rs-fMRI map intensity + MRI rating (with balanced accuracy = 0.71). Only the latter model (rs-fMRI map intensity + MRI rating + clinical variables) resulted to be significantly different from the model with only clinical variables (with balanced accuracy = 0.64 and AUC = 0.71) in terms of classification results (McNemar test, p < 0.05).
Statistically, accuracy and balanced accuracy had near value because the sensitivity and specificity of AUC were similar. Considering the significant variables included by LASSO method there were no significant differences between right and left hemispheres ( Table 5 and Supplementary Table S6).

Differences in Etiology
Considering possible differences in etiology, we again conducted logistic regressions analyses after applying the LASSO method (Table 6, Supplementary Tables S7 and S8). In traumatic patients, the LASSO method did not retain any variable. In vascular patients, the balanced accuracy was 0.74 for both rs-fMRI indices. Combining rs-fMRI indices with MRI rating, the results did not change; adding the clinical variables, the accuracy did not improve and only the disease duration was significant as a variable. In post-anoxic patients, the balanced accuracy was 0.63 for rs-fMRI rating and 0.90 for rs-fMRI map intensity. Combining the two rs-fMRI indices with MRI rating, the balanced accuracy became 0.63 and 0.83, respectively. Adding the clinical variables, the balanced accuracy did not change considerably and only the disease duration was significant as a clinical variable (see Supplementary Table S9). Only for this class, LASSO included more variables of the left than the right hemisphere.
Overall, considering the three etiological classes, all networks contributed similarly to accuracy.
Considering the number of networks identified, a positive relation was observed with the CRS-R total score for the total of patients (ρ = 0.26, p < 0.01), indicating that as the CRS-R total score increases, the number of networks also increases ( Table 7); this relation was stronger for the anoxic class (ρ = 0.39, p < 0.05). Generally, the two visual networks were the most frequently recognizable components. In particular, for traumatic patients, the number of networks identified generally increased as the diagnostic group passed from VS/UWS to MCS to SD ( Table 7). The visual networks better distinguished VS/UWS from MCS compared to the other components. For vascular patients the recognition of networks did not increase with diagnosis, rather the effect was opposite, with a higher number of networks for VS/UWS than MCS. Only the AUD network was more frequently observed for MCS than VS/UWS. This might be due to the lesion site, often involving the territory of the middle cerebral artery in vascular injury, with a more severe impact for VS/UWS than MCS. For anoxic patients, the recognition of networks clearly increased as the diagnosis improved; the visual networks were the less recognizable component for VS/UWS (13%), but the more recognizable for MCS (80%).
Finally, the Spearman correlation between the number of networks observed for VS/UWS and MCS (without SD) and CRS-R total score was not significant, nor was the one with CRS-R-Modified score.
Furthermore, there was a negative relation between the number of networks identified and disease duration for the total of patients (ρ = −0.19, p < 0.05), indicating that as disease duration increases, the number of networks decreases (Table 7).
Clinical data of all patients are reported in Supplementary Table S10.

Discussion
In 109 patients with DoC and emerged from DoC we investigated the integrity of the four rs-fMRI networks more related to primary functions, namely the SM, LVIS, MVIS and AUD components, and analyzed data in relation to structural MRI and clinical information (etiology, disease duration and age), as occurs in clinical practice.

Assessment of Single Networks
Our data show that: (i) VS/UWS patients present a lower number of networks compared to MCS patients; (ii) functional connectivity is useful to differentiate VS/UWS from MCS; and (iii) the MVIS and LVIS networks are correlated with the clinical status measured with CRS-R.
In particular, VS/UWS patients were most frequently associated with 0 networks and MCS patients to 2 networks, although both classes ranged from 0 to 4 components. Only SD patients never had any networks and were most frequently associated with 4 components.
The diagnostic accuracy provided by single networks is generally limited when the rs-fMRI data is taken alone (on average balanced accuracy = 0.57 and AUC = 0.55). To represent the clinical reality where rs-fMRI is typically assessed with structural imaging, we included a model with rs-fMRI+MRI data from each network, which however showed low accuracy values (on average a balanced accuracy of 0.58 and AUC = 0.67). As these are lower-order rs-fMRI networks corresponding to primary functions, an inferior diagnostic power is not unexpected compared to higher-order networks such as DMN.
Results are lower than those reported by Demertzi et al. (2015) [25] which demonstrated that visual, auditory and sensorimotor networks assessed individually had a discriminative capacity greater than 80%; however, they studied a smaller sample of patients in a subacute phase and used different statistical methods. To date, no other studies have investigated these networks with respect to the VS/UWS vs. MCS classification and differentiation.
The two rs-fMRI measures used here are rs-fMRI rating and rs-fMRI map intensity: the first one is based on a qualitative, time-consuming assessment indicating the degree of network recognition, probably leading to more interpretable results; the second one results from an automated analysis procedure indicating the mean intensity of the network. Interestingly, both measures yielded similar results, which demonstrates that implementing rs-fMRI into the clinical routine based on expert assessment can be feasible [24], although automatic rating can be quicker especially when the sample size is large.
A potentially relevant question for the clinical assessment refers to the correspondence between the rs-fMRI network and the CRS-R subscale and is whether the presence of a visual, auditory or sensorimotor map can be clinically useful for evaluating patients' perceptual functioning. In our sample, only the visual subscale could be correlated with the rs-fMRI data, as scores of the motor and auditory CRS-R subscales were centered on a few values. Results showed that the highest the visual score was, the more recognizable resulted the rs-fMRI visual map, indicating that the type of response and the level of behavioral complexity in the visual domain were correlated with the integrity of LVIS and MVIS networks. The voxel-wise analysis also showed a significant correlation between the integrity of the visual networks and the CRS-R visual score: as visual scores > 2 are associated with MCS, a well-preserved rs-fMRI visual map might be associated with a high level of consciousness. This, however, does not imply that the presence of the visual network is a prerequisite for, or a direct sign of consciousness; these correlative results only suggest that the visual networks can provide additional information on the visual functioning in DoC and might be functionally relevant for the state of consciousness of the patient.
More generally, the correspondence between CRS-R subscales and network recognizability allowed us to assess those cases where the clinical score was zero, i.e., no behavioral response, but the rs-fMRI map was present. Obviously, they can be cases of possible overestimation, where a network is erroneously identified, despite a low clinical score; alternatively, they can represent a concrete demonstration of how functional imaging measurements can reveal the presence of residual brain functioning, despite scarce auditory/sensorimotor/visual/cognitive performance. For some cases reported in the literature, the presence of a rs-fMRI network map indicated residual brain function, not detectable at the patient's bed and associated with clinical recovery [24,52]. In our study, 7 cases had a clinical score = 0 but showed a rs-fMRI network that was clearly recognizable, which was confirmed by EPs and FDG-PET. This finding confirms that the presence of clearly identifiable rs-fMRI networks in patients with DoC may warn clinicians about possible clinical assessment biases and prompt a more careful evaluation of the patient's condition. Indeed, in some patients, motor impairments can cause difficulties in clinical assessment [53]. In this context, the rs-fMRI technique is a valuable supplementary assessment [8] that can be added to multimodal DoC evaluations, providing information on residual cortical functions.

Assessment of Four Networks
Our data show a significant relationship between the number of networks identified and the clinical status measured with the CRS-R total score (ρ = 0.26, p < 0.01). Combining data from the four networks, rs-fMRI map intensity achieved the highest accuracy scores (balanced accuracy = 0.67 and AUC = 0.77) which were even higher when the model included also MRI, yielding very good classification results (balanced accuracy = 0.71 and AUC = 0.80) [51]. Importantly, the model with rs-fMRI map intensity + MRI + clinical variables was significantly better than the one considering only clinical variables. Generally, all networks contributed to the results. We, therefore, advise that rs-fMRI and structural imaging should be considered together to attain the best diagnostic accuracy.

Effect of Etiology and Disease Duration
Two issues, scarcely investigated yet, are the effects of etiology and disease duration in patients with DoC. The ability to classify VS/UWS vs. MCS was assessed separately for the traumatic, vascular and anoxic etiologies. Balanced accuracy and AUC were higher for anoxic and lower for traumatic etiology, while for the vascular class, more networks were apparently observed for VS/UWS than MCS. These results confirm that the VS/UWS vs. MCS distinction is more difficult for the vascular (ischemic and hemorrhagic) etiology, where the damage to cortex and white matter is typically more nuanced and heterogeneous, and clusters of correlated activity are more difficult to judge as neural or artifactual [24,39]. For the anoxic class, the rs-fMRI map intensity had a balanced accuracy of 0.90 and AUC = 0.98, indicating certain ease of classification for this etiology despite the high clinical complexity [31]; the result was confirmed by a significant relationship between the number of networks observed and CRS-R total score (ρ = 0.39, p < 0.01). For traumatic patients, there was no discrimination between VS/UWS vs. MCS patients, probably due to the location of the damage which, on average, made accuracy more difficult. Overall, results warn about the impact of etiology on the analysis of rs-fMRI networks where vascular patients are more difficult to classify and suggest subdividing patients by etiology if the patient sample allows it.
Among the four networks, the visual components were the most frequently recognizable in particular in traumatic and vascular patients. This is in agreement with recent work on sedated and unsedated patients with DoC showing that among seven rs-fMRI networks the visual network was the most similar among controls and unsedated patients [14].
Disease duration had a weak but significant impact on the global number of networks observed (ρ = −0.19, p < 0.05), where, as duration increased, the number of networks decreased; when diagnostic accuracy was considered in relation to the different etiologies (Table 4), LASSO method considered disease duration the only significant clinical variable to combine with imaging data (Table S2 Supplementary Material). However, disease duration had no effect when networks were analyzed individually.
Finally, we would like to discuss some methodological issues of rs-fMRI analysis in DoC. Several important choices can affect the recognizability of networks in patients. These include the global signal regression, the motion parameter regression, and the number of components in ICA. Regarding the motion parameter regression, recent findings suggest that it is not necessary to have a high number of regression parameters (e.g., 24 parameters or more) to minimize motion artifacts [54]. In addition, censoring strategies, such as despike and spike regression, significantly improved network identifiability [54,55], which is essential in DoC. Regarding the global signal regression, there is no consensus on its inclusion in rs-fMRI analyses, even in DoC. Some studies have included it [27,47], whereas others did not [25,56,57]. Moreover, global signal regression may also increase anticorrelations [58] and, when motion artifacts are present, it can introduce connectivity patterns in which connectivity increases when nodes are close to each other, and decreases as the distance between nodes decreases [54,55]. For these reasons, our analysis included the regression of outlier volumes (spike regression) and of a limited number of motion parameters (3 rotation + 3 translation + WM signal + CSF signal) without global signal regression.
With regard to the selection of ICA components, the choice of using automatic methods or a fixed number of components can further impact the identification of networks. An automated method may result in a large number of components, leading to the splitting of networks, or conversely, a low number of components could lead to an underestimation of the number of networks. The choice of a fixed number of components (e.g., 20-30) is a good compromise to identify all the neural components of interest as well as to identify the known artifactual components in these patients. In addition, this type of patient has severe neuroanatomical alterations and it is often difficult to determine whether a component is neural or vascular. To better distinguish the neural components from noise and artifact, we first performed an automatic identification and then we visually verified the components, according to Griffanti's flowchart [34]. Finally, the choice of a fixed number of components (e.g., [20][21][22][23][24][25][26][27][28][29][30] is in line with other rs-fMRI studies on DoC [24,29,43,59].

Limitations and Future Directions
In this study, we did not consider follow-up clinical data, as it was only sparingly available. Our group of patients had non-heterogeneous scores in the CRS-R subscales, except for the visual function and this represents a limitation for rs-fMRI data correlations. In the post-anoxic etiology, there was an unbalanced proportion of VS vs. MCS (30 vs. 5), which is however typically observed in the clinical realm. For the diagnostic accuracy of the four networks subdivided by etiology, traumatic patients were associated with a low diagnostic accuracy. Future work can group traumatic patients with similar damage locations. For the rs-fMRI quantitative analysis we calculated the mean intensity of the network nodes. In the future other dimensions can be explored, for example integrating the spatial information of the connectivity map with the mean intensity. The networks included in this study are closely related to specific functions. Future studies should shed light on the comparison between accuracy provided by lower-level vs. higher-level networks.

Conclusions
The four rs-fMRI networks more related to primary functions may provide additional information on the visual, auditory and sensory-motor pathways, in some cases comparable to EP and FDG-PET techniques and in absence of clinical responses. The MVIS and LVIS networks are correlated with the clinical status, as the CRS-R subscale scores of our patients were more distributed, unlike the auditory and motor subscales.
The four rs-fMRI networks have limited diagnostic accuracy when assessed individually but, when considered altogether with MRI, they increase to a very good level (AUC = 0.80) and the model rs-fMRI map intensity+MRI rating provides additional significant information to the clinical data.
Generally, rs-fMRI analyses based on a qualitative rating or a quantitative map intensity provide convergent results. Among the clinical variables, etiology can affect diagnostic accuracy and disease duration is a significant factor in the models examined.
In conclusion, rs-fMRI is a non-invasive technique, indicated when an MRI exam is performed [8] and our data confirm that it can be clinically useful in the assessment of patients with DoC.
Supplementary Materials: The following supporting information can be downloaded at: https: //www.mdpi.com/article/10.3390/brainsci12030355/s1, Table S1: AUC, confidence interval of AUC, sensitivity, specificity, and patients wrongly classified as VS/UWS and MCS by etiology as traumatic/vascular/anoxic of single networks reported in Table 2 using LOOCV; Table S2: AUC, confidence interval of AUC, sensitivity, specificity for each single network of five models using 10-fold CV; Table S3: Variables included by LASSO method for diagnostic accuracy of single networks reported in Table 2; Table S4: AUC, confidence interval of AUC, sensitivity, specificity, and patients wrongly classified as VS/UWS and MCS by etiology as traumatic/vascular/anoxic of 4 networks reported in Table 5 using LOOCV; Table S5: AUC, confidence interval of AUC, sensitivity, specificity for the 4 models using 10-fold CV; Table S6: Variables included by LASSO method for diagnostic accuracy of the 4 networks reported in Table 5; Table S7: AUC, confidence interval of AUC, sensitivity, specificity, and patients wrongly classified as VS/UWS reported in Table 6 using LOOCV; Table S8: AUC, confidence interval of AUC, sensitivity, specificity for the 4 models by etiology using 10-fold CV; Table S9: Variables included by LASSO method for diagnostic accuracy of the 4 networks considering the different etiologies reported in Table 6; Table S10: Clinical data are reported for all patients.

Institutional Review Board Statement:
The study was conducted in accordance with the Declaration of Helsinki, and approved by the Ethics Committee of Fondazione IRCCS Istituto Neurologico Carlo Besta (Resolution N IX/000407 of 05/08/2010).

Informed Consent Statement:
Informed consent was obtained from the legally authorized representative of the patients prior to inclusion in the study.

Data Availability Statement:
The data presented in this study are available on request from the corresponding author.

Conflicts of Interest:
The authors declare no conflict of interest.