Alteration of brain network topology in HIV-associated neurocognitive disorder: A novel functional connectivity perspective.

HIV is capable of invading the brain soon after seroconversion. This ultimately can lead to deficits in multiple cognitive domains commonly referred to as HIV-associated neurocognitive disorders (HAND). Clinical diagnosis of such deficits requires detailed neuropsychological assessment but clinical signs may be difficult to detect during asymptomatic injury of the central nervous system (CNS). Therefore neuroimaging biomarkers are of particular interest in HAND. In this study, we constructed brain connectivity profiles of 40 subjects (20 HIV positive subjects and 20 age-matched seronegative controls) using two different methods: a non-linear mutual connectivity analysis approach and a conventional method based on Pearson's correlation. These profiles were then summarized using graph-theoretic methods characterizing their topological network properties. Standard clinical and laboratory assessments were performed and a battery of neuropsychological (NP) tests was administered for all participating subjects. Based on NP testing, 14 of the seropositive subjects exhibited mild neurologic impairment. Subsequently, we analyzed associations between the network derived measures and neuropsychological assessment scores as well as common clinical laboratory plasma markers (CD4 cell count, HIV RNA) after adjusting for age and gender. Mutual connectivity analysis derived graph-theoretic measures, Modularity and Small Worldness, were significantly (p < 0.05, FDR adjusted) associated with the Executive as well as Overall z-score of NP performance. In contrast, network measures derived from conventional correlation-based connectivity did not yield any significant results. Thus, changes in connectivity can be captured using advanced time-series analysis techniques. The demonstrated associations between imaging-derived graph-theoretic properties of brain networks with neuropsychological performance, provides opportunities to further investigate the evolution of HAND in larger, longitudinal studies. Our analysis approach, involving non-linear time-series analysis in conjunction with graph theory, is promising and it may prove to be useful not only in HAND but also in other neurodegenerative disorders.

Introduction

It is estimated that > 50% of the individuals affected with HIV are susceptible to developing HIV-associated neurocognitive disorder (HAND) (Ellis et al., 2007, Heaton et al., 2015, Masters and Ances, 2014). The virus is capable of crossing the blood brain barrier at an early stage during infection, establishing an inflammatory milieu that causes damage to synaptodendritic connections resulting in neuronal dysfunction (Ellis et al., 2007). Such changes at the neuronal level can occur much earlier than the cognitive symptoms in an individual, (Cysique and Brew, 2011, Holt et al., 2012, Wang et al., 2011). With an increased life-expectancy owing to therapeutic advancements in treatment of HIV positive (HIV +) individuals, neurocognitive deficits affecting multiple domains eventually develop in HIV infected individuals (Antinori et al., 2007b, Thomas et al., 2013). Thus, given the high likelihood of developing HAND in the long-term, there is an increasing demand for biomarkers that have the potential to improve detection of injury of the central nervous system (CNS).

In this context, neuroimaging can serve as a potential tool for providing quantitative insights into the progression of HAND. Imaging techniques such as magnetic resonance spectroscopy, volumetrics, diffusion tensor imaging and PET can be used to study this disease and develop a better understanding of its progression (Holt et al., 2012). Using quantitative MRI techniques, recent studies (Jahanshad et al., 2012, Masters and Ances, 2014, Zhu et al., 2013), have demonstrated thinning of cortical gray matter and specific disruption to white matter integrity within the brain, in subjects with HAND. Further, the use of BOLD functional MRI signal acquired in the absence of external stimulation can give novel insights into the complex organization of the brain (Lee et al., 2013). Few recent studies (Thomas et al., 2013, Wang et al., 2011) have shown alterations in the interactions of specific subsystems such as the lateral occipital cortex (LOC), default mode network (DMN) and salience (SAL) networks in subjects with HIV infection. Recently, (Ortega et al., 2015) have shown that connectivity between cortico-striatal networks is affected in HIV + individuals, however these changes were not correlated with CD4 T-cell counts, plasma viral loads or neuropsychological assessment scores. Changes to characteristics of the brain network during the course of HIV infection has also been explored using magnetoencephalography (Wilson et al., 2015). Although efforts have been made employing different techniques, newer methods are required to develop a better understanding of the pathophysiology of HAND (Saylor et al., 2016).

The use of graph theory in conjunction with fMRI analysis, has offered interesting insights into the organizational principles of the complex integrative brain network (Bullmore and Sporns, 2009, Bullmore and Bassett, 2011, Fornito et al., 2013). Application of such analysis techniques has shed light on characteristic changes occurring in different neurological disorders (Stam, 2014). For example, it has been shown using graph theoretic approaches that there is a subtle randomness occurring in subjects which schizophrenia (Lynall et al., 2010), loss of ‘small-world’ properties in Alzheimer's disease (Sanz-Arigita et al., 2010, Stam et al., 2007), significant changes in network modularity in Parkinson's disease (Baggio et al., 2014) etc. However, not many studies (Thomas et al., 2015) have investigated HAND from the perspective of graph theory.

Moreover, almost all earlier studies investigating functional connectivity changes in HIV using resting state fMRI have used Pearson's correlation based approaches on low-pass or band-pass filtered fMRI time series to quantify their functional connectivity. The use of correlation, though simpler in application, tends to ignore the non-linear characteristics of the brain (Daunizeau et al., 2012, Zhang et al., 2008). We have previously proposed Mutual Connectivity Analysis (MCA) (Wismüller et al., 2016, Wismüller et al., 2014) as an alternative to correlation to capture the intrinsic dynamics of the brain networks. In contrast to correlation, this approach is capable of quantifying non-linear interdependence between brain time-series. We hypothesize that such a technique could potentially reveal additional information about the brain networks which may be critical to predict and characterize the progression of neurological disorders.

This study presents an analysis of HAND using advanced fMRI time-series analysis techniques. We aimed to test whether MCA could uncover differences in the topological properties of the brain in HIV infected individuals compared to healthy subjects and also analyze their ability to supplement information provided by the Neuropsychological (NP) testing derived z-scores.

Materials and methods

Study participants

Subjects participating in this study were recruited at the University of Rochester Medical Center, as part of a NIH funded study (R01-DA-034977). All participants provided written informed consent prior to participation as per protocol approved by the institutional IRB. Exclusion criteria for the study included any pre-existing psychiatric or neurological disorder, active brain infection (except HIV-1), presence of brain neoplasm or space-occupying lesion, chronic seizures or head injury, active substance abuse and/or other criteria which made subject ineligible for MR imaging. The serostatus of HIV + individuals was confirmed by detection of plasma HIV RNA or by documented positive HIV enzyme-linked immunoassay. The final cohort consisted of a total of 20 HIV + subjects (mean age 42 ± 10 years; 9 females) and 20 age-matched control subjects (mean age 45 ± 16 years; 5 females), who were selected based on the absence of confounding neurological disorders and a normal neurological exam. All participating subjects underwent an extensive standard battery of neuropsychological (NP) tests covering 7 cognitive domains: verbal fluency, executive function, speed of information processing, attention/working memory, learning, and memory motor function (perceptual-motor speed). The specific tests performed are detailed in Supplementary Table 1. The raw NP test scores were converted to summary demographically adjusted z-scores of cognitive performance across multiple cognitive domains. The domain z-scores were used for the assessment of HAND (Antinori et al., 2007a) and studying associations with imaging derived parameters. The use of the overall Z-score is often used as the primary cognitive outcome compared to single cognitive domains because of its sensitivity (Schifitto et al., 2001, Schifitto et al., 2007).

The exact steps followed in our analysis are summarized in Fig. 1.

Fig. 1

Schematic representation of the different steps in the study.

Imaging protocols

All subjects underwent MRI scanning, using a 3.0 Tesla Siemens Magnetom TrioTim scanner equipped with a 32 channel phased array head-coil, at the Rochester Center for Brain Imaging (Rochester, NY, US). A T1-weighted magnetization-prepared rapid gradient acquisition echo sequence (MPRAGE) was performed with the following parameters: repetition time = 2530 ms, echo time = 3.44 ms, isotropic voxel size 1 mm, flip angle = 7°, acquisition time ~ 3 min. This sequence was used mainly for co-registration of the functional sequence to the standard MNI152 template.

Functional scanning was performed using echo planar imaging (EPI) with the following parameters: repetition time = 1650 ms, echo time = 60 ms, 25 contiguous slices with a 96 × 96 acquisition matrix and 2.5 mm resolution, flip angle = 84°. A total of 250 volumes were acquired during a scan of ~ 6 min 40 s. The subjects were instructed to lie down still with their eyes closed and not think about anything particular.

Preprocessing

Prior to computation of connectivity, the fMRI data used in this study was preprocessed using standard methodology. Briefly, the following operations were performed, for each dataset, the first ten (of 250) volumes of functional magnetic resonance images were removed to analyze only those in which steady-state imaging had been reached. Next, motion correction, brain extraction and correction for slice timing acquisition were performed. Additional nuisance regression, to remove variations due to head motion and physiological processes was carried out. This model included linear and quadratic trends, signals from white matter and cerebrospinal fluid using CompCor (Behzadi et al., 2007), and the Friston-24 motion parameters (Friston et al., 1996). Each individual's dataset was registered to the 2 mm MNI standard space using a 12-parameter affine transformation. All preprocessing steps were carried out using the C-PAC software (Sikka et al., 2014) and its corresponding dependencies in FSL (Smith et al., 2004).

We applied the commonly used Automated Anatomic Labeling (AAL) template (Tzourio-Mazoyer et al., 2002), for defining ROIs in the cortical and subcortical regions of the brain. A total of 90 regions (45 belonging to each hemisphere excluding the brain stem and cerebellar regions) were defined. These regions form nodes of our system, and representative time-series for each ROI were computed by averaging the time-series of all voxels within it.

Functional connectivity

Mutual Connectivity Analysis can be used to detect non-linear functional relationships between pairs of time-series. The pair-wise affinity between regional time-series describes the degree of their dynamic coupling as a measure of their cross-prediction performance. For two time series, x and y of length l, we first computed the influence θ of x on y by constructing time-delayed vectors x and y with embedding dimensions m and n, respectively, i.e. x = ((x(t), …, x(t + m-1))T and y = ((y(τ), …, y(τ + n-1))T, where τ = t + m and t ∈ {1,2, …, l-m-n + 1}. These were paired such that each y was associated with a vector x that temporally precedes it. For example, when embedding dimensions are chosen as m = 5 and n = 1, then x1 = (x(1), x(2), …, x(5))T and y1 = (y(6)) form a pair. We used a generalized radial basis function (GRBF) neural network to approximate a non-linear mapping Φ. We adopted a random subsampling cross-validation strategy with 70% training (Tr) and 30% test (Te) data over 10 iterations. We used the set of pairs {(xTr,yTr)} as training data for computing the non-linear mapping Φ, which was then used on an independent the test dataset {xtTe} for prediction. The correlation between the predicted and the true target vectors yTe was averaged over the 10 iterations as a measure θ of the influence of x on y based on non-linear predictability. This procedure was repeated for all pairs of time-series x and y in the dataset. More detailed information on the computations can be found in the supplementary material and our previous work (Wismüller et al., 2016, Wismüller et al., 2014).

In the current dataset, we have N = 90 time-series. The directional interactions for all N2 pairs were computed and stored in an affinity matrix θ. The directional matrix was symmetrized (A = (θ + θT)/2) to produce an undirected network for further graph-theoretic analysis. The self-connections for the individual nodes were not considered for further analysis.

Additionally, to study graph-theoretic properties of networks using a conventional method, we low-pass filtered the 90 regional time-series at 0.1 Hz in line with numerous prior studies suggesting that endogenous fMRI dynamics of neuronal origin are captured below this frequency (Lynall et al., 2010, Sanz-Arigita et al., 2010). These time-series were then used to construct a network connectivity profile P based on Pearson's correlation for all regional pairs.

Graph-theoretic analysis

The topological properties of the brain network derived from both methods our MCA-GRBF and the conventional Pearson's correlation approach were defined on the basis of a graph defined as follows.where e refers to an existing connection edge from node i to node j in the graph, if the connection strength between the nodes exceeds a particular threshold, φ. The choice of φ for this study was based on connection density which is defined as the ratio of the number of edges to the total number of possible edges in the graph. This ensures an equal number of edges the graphs obtained for each subject, increasing the reliability of inter-subject comparisons. The appropriate threshold was chosen based on the specific criteria as outlined in previous work (Lynall et al., 2010), i.e. the graphs should be fully connected and have non-random organization (Humphries et al., 2006, Kaiser and Hilgetag, 2006). Using these conditions, we obtained range of densities between 0.40 and 0.58. As a similar trend for intergroup differences was observed over this range of thresholds, we report the results obtained chose a threshold preserving 40% of the connections for further analyses in line with prior work (Baggio et al., 2014).

To characterize the overall topological organization of the brain network matrices, the following commonly used graph-theoretic metrics were computed. These are discussed here briefly. Extensive descriptions are available elsewhere e.g. (Bullmore and Sporns, 2009, Rubinov and Sporns, 2010):

Degree represents the number of links an individual node is connected to. The variance of degree across all nodes indicates the homogeneity of the overall network. Network integration, a property that represents the ability to combine the information from distributed areas, is well captured by measuring the average path length between the various regions. In this study, Efficiency, which is proportional to the inverse of characteristic path length, was computed as a measure of network integration. Network segregation metrics are used to detect the presence of specialized regions (or clusters) within a network indicating segregated functional dependencies in the brain. This is represented by the Clustering Coefficient which is the fraction of neighboring nodes being also nodes to each other. Another useful measure, Modularity, is used to describe the degree to which a network can be divided into clearly delineated and non-overlapping sub-groups. The brain is known to simultaneously reconcile the opposing demands of functional integration and segregation, which can be captured using the measure, Small-Worldness. The resilience of a brain network represents its ability to withstand alterations of its overall activity which may be induced by neuropathological damage. The Assortativity coefficient is commonly used as a measure of resilience. Betweenness Centrality of each node in the network to explore the role of certain regions to function as central ‘hubs’ in a network (Caeyenberghs et al., 2012) and analyze, if these change with respect of HAND progression.

Statistical analysis

Comparisons between two independent groups were conducted by two-group Welch's unequal variances t-test. Multiple regression analyses were used to quantify the linear associations between continuous variables derived from fMRI analysis, neuropsychological tests as well as other clinical variables while controlling for age and gender effects. In these analyses, linear coefficients were estimated by the maximum likelihood criterion, and the statistical significance was assessed by regression t-test. A p-value p < 0.05 was considered statistically significant for a single hypothesis testing problem. For inferential problems that involved multiple hypotheses, Benjamini–Hochberg multiple testing procedure (Benjamini and Hochberg, 1995) was used to control the false discovery rate (FDR) at < 0.05 level. Briefly, this procedure strives to control the mean ratio of false discoveries among all rejected hypotheses at a given threshold, e.g., 0.05. Consequently, we can expect that on average, < 5% of all significant discoveries (those hypotheses with adjusted p-values < 0.05) are false discoveries. All statistical analyses were performed in R 3.3.0 (R Foundation for Statistical Computing, Vienna, Austria).

Results

Demographics and descriptive statistics

Table 1 summarizes the demographic and clinical data of the subjects in this study. There were no significant differences in age between the two subject groups (p > 0.05), though HIV − control group contained a larger number of female subjects. In seropositive individuals, absolute CD4 counts ranged from 74 to 1730/mm3, whereas plasma viral load ranged < 40/mL (undetectable) to 62,000 copies/mL. With the exception of one individual, all HIV + subjects were on a stable cART regimen. Of these, four of the subjects showed viral loads > 20,000 copies/mL at the time of the visit. Viral loads were undetectable (< 40) in 11 of the subjects. The duration of HIV infection within our cohort, calculated based on the time of initial HIV infection diagnosis, ranged from a few months to 16 years. This however may not reflect the true duration of infection.

Table 1

Demographic details and clinical characteristics for the population used in this study. p-values indicate differences between the HIV + and control subjects.

	HIV −				HIV +				p-Value
	Mean	(SD)	Range		Mean	(SD)	Range
Number of patients	20				20
Age – in years	41	(10)	21	60	42	(15)	23	71	> 0.05
Gender (female/male)	9/11				5/15				-NA-
Nadir CD4 (cells/mm3)	-NA-				314.5	(212)	12	710	-NA-
CD4 (cells/mm3)	-NA-				703	(465)	74	1730	-NA-
VL (log10 scale)	-NA-				1.76	(1.79)	0	4.79	-NA-
HIV - in years	-NA-				11	(9)	0.08	26	-NA-
NP Z-scores
Attention	0.432	(0.8)	− 1.20	1.60	− 0.514	(1.021)	− 2.851	1.274	0.002
Executive	0.284	(0.982)	− 1.22	2.57	− 0.36	(0.953)	− 1.884	1.636	0.079
Learning	0.363	(0.915)	− 1.39	1.81	− 0.323	(0.878)	− 1.782	1.361	0.015
Memory	0.345	(1.04)	− 1.79	1.86	− 0.254	(0.760)	− 1.689	1.263	0.036
Motor	0.536	(0.664)	− 0.49	2.09	− 0.54	(0.957)	− 1.982	1.615	0.001
Speed of information processing	0.366	(0.769)	− 1.50	1.81	− 0.42	(1.066)	− 2.328	1.269	0.041
Overall	2.327	(2.842)	− 3.27	6.97	− 2.41	(3.711)	− 8.407	6.435	< 10− 4
HAND classification (%)
WNL	-NA-				6	(30%)			-NA-
ANI	-NA-				12	(60%)			-NA-
MND	-NA-				2	(10%)			-NA-
HAD	-NA-				0	(0%)			-NA-

NP testing and HIV status association

Cognitive performance of HIV infected subjects as quantified by their NP z-scores (Table 1), are invariably worse than the healthy controls in all cognitive domains. HIV infection has negative effects on all specific domain z-scores, with the age and gender effects being controlled for. Regression t-tests revealed that these differences are significant in all but the executive domain (p = 0.079), while the most significant impact is observed on the overall z-score (p < 10− 4). Significant confounding effects of age and gender are seen for executive, spatial, attention and motor domains. These results are summarized in Supplementary Table 2.

The recommendations as defined in (Antinori et al., 2007a) were used to define HAND status within HIV + subjects into four categories in our cohort – ANI (asymptomatic neurological impairment), MND (mild neurological disorder), WNL (within normal limits), and HAD (HIV associated dementia). No subjects met the criteria for the most severe neurological damage i.e. HAD. Details of the number of subjects in each class are provided in Table 1.

Connectivity profiles

The mean connectivity profiles obtained using the MCA-GRBF approach as well as correlation for the subject cohorts are shown in Fig. 2. Although an overall similarity of patterns is observed we note certain differences in some regions. Further graph theoretic-analyses were carried out to quantify such differences.

Fig. 2

Mean connectivity profile for the control subjects (left) and patients (right). Gross labels for the regions have been provided here, for clarity. Exact region names (in order) are listed in Supplementary Table 2. The figure shows that matrices obtained using MCA-GRBF and correlation seems to capture different information with regard to the connectivity between the different brain regions. Subtle differences (inferior parietal to frontal lobe connectivity (p < 0.05), connectivity within the cingulate regions (p < 0.1); MCA-GRBF) are seen when comparing the mean matrix of HIV + subjects and controls, and these are further captured by using graph theoretic analysis in the subsequent section(s).

Relationship between imaging metrics and NP testing

We tested the relationships between imaging derived measures and cognitive performance z-scores based on NP testing (which is the current clinical standard for diagnosis of HAND), an initial exploratory ordinary linear regression analysis showed that there are significant marginal associations between the two sets of parameters (e.g. overall z-score with Modularity and Small-Worldness). A few of these marginal associations are visualized in Supplementary Fig. 1. Further investigation based on multiple linear regressions showed significant associations between MCA-GRBF derived imaging parameters and NP tests (detailed in Table 2).

Table 2

p-Values (uncorrected and with FDR correction) for the multiple regression models of image derived variables with NP testing z-scores. Significant associations (p < 0.05) are highlighted. (CC-Clustering co-efficient, DV-Degree-Variance, SW-Small-Worldness).

	NP-zscore	ATTENTION		EXECUTIVE		LEARNING		MEMORY		MOTOR		SPEED		OVERALL
	Measure	p	p.fdr	p	p.fdr	p	p.fdr	p	p.fdr	p	p.fdr	p	p.fdr	p	p.fdr
Cross correlation	Assortativity	0.701	0.789	0.668	0.785	0.305	0.549	0.085	0.284	0.050	0.284	0.403	0.659	0.477	0.745
	CC	0.487	0.745	0.056	0.284	0.166	0.391	0.606	0.759	0.253	0.471	0.086	0.284	0.080	0.284
	DV	0.180	0.391	0.005	0.150	0.084	0.284	0.331	0.577	0.098	0.284	0.027	0.284	0.010	0.171
	Efficiency	0.100	0.284	0.063	0.284	0.598	0.759	0.888	0.905	0.851	0.884	0.549	0.759	0.228	0.439
	Modularity	0.086	0.284	0.015	0.196	0.066	0.284	0.117	0.305	0.062	0.284	0.086	0.284	0.006	0.150
	SW	0.619	0.759	0.181	0.391	0.595	0.759	0.921	0.921	0.175	0.391	0.119	0.305	0.208	0.433
MCA - GRBF	Assortativity	0.017	0.091	0.193	0.326	0.584	0.671	0.681	0.750	0.025	0.091	0.336	0.432	0.065	0.148
	CC	0.230	0.354	0.025	0.091	0.303	0.408	0.580	0.671	0.157	0.292	0.049	0.127	0.042	0.117
	DV	0.033	0.105	0.001	0.019	0.324	0.427	0.779	0.789	0.044	0.117	0.022	0.091	0.008	0.070
	Efficiency	0.357	0.448	0.036	0.108	0.789	0.789	0.778	0.789	0.691	0.750	0.264	0.376	0.283	0.392
	Modularity	0.016	0.091	0.000	0.019	0.211	0.336	0.560	0.671	0.025	0.091	0.029	0.097	0.003	0.038
	SW	0.020	0.091	0.001	0.019	0.260	0.376	0.752	0.789	0.015	0.091	0.012	0.091	0.004	0.038

Overall z-score is a significant predictor for measures based on matrices obtained using the MCA-GRBF (Modularity - p = 0.0377, Small-Worldness - p = 0.0384). Additionally, highly significant effects are seen when using executive z-score as a co-variate for Modularity, Small-Worldness and Degree-Variance. The traditional correlation based metrics, however, do not produce any significant associations with clinical z-score derived from of NP testing.

Additionally, we performed a similar analysis, segregated by the two different cohorts. As it could be expected, by separating the population into two cohorts, the spread of cognitive scores to correlate with the graph metrics is decreased, especially in the control group who tend to have more homogenous cognitive performance than the HIV infected group. This was reflected in fewer significant associations (p < 0.05, uncorrected) between the metrics within control group than in the HIV group (Supplementary Table 3).

Effect of HIV status on imaging metrics

We additionally studied the effect of HIV status on our imaging derived characteristics. We performed multiple regression analyses to study the linear associations between graph-theoretic measures (as response variables) and HIV infection status (the main covariate). Regression t-tests were used to determine the significance of the HIV effects. Mean (SD) of graph-theoretic measures in both cohorts and the p-values of the HIV effects are summarized in Table 3.

Table 3

Image derived parameters of subjects with and without HIV infection. p-values (with FDR correction) indicate significance of association of the parameters with HIV status.

	HIV − mean (SD)		HIV + Mean (SD)		p-Value (uncorrected)	p-Value (FDR corrected)
MCA-GRBF derived measures
Assortativity	− 0.041	(0.116)	0.011	(0.114)	0.097	0.201
Clustering Coefficient	0.857	(0.047)	0.841	(0.036)	0.212	0.336
Degree-Variance	473.302	(90.474)	417.492	(110.242)	0.068	0.148
Efficiency	0.670	(0.033)	0.686	(0.038)	0.157	0.292
Modularity	0.131	(0.036)	0.155	(0.05)	0.065	0.148
Small-Worldness	1.124	(0.075)	1.173	(0.117)	0.114	0.228
Correlation derived measures
Assortativity	0.113	(0.067)	0.144	(0.084)	0.224	0.439
Clustering Coefficient	0.813	(0.04)	0.804	(0.042)	0.708	0.789
Degree-Variance	331.081	(77.811)	298.066	(84.23)	0.353	0.597
Efficiency	0.719	(0.023)	0.724	(0.01)	0.562	0.759
Modularity	0.173	(0.025)	0.191	(0.034)	0.100	0.284
Small-Worldness	1.213	(0.068)	1.235	(0.078)	0.516	0.753

The functional networks derived using MCA-GRBF expressed some key small-world characteristics consistently across both cohorts (e.g.: Small Worldness– control group: mean = 1.124, SD = 0.075; HIV + group: mean = 1.173, SD = 0.117), similar to what has been noted in earlier studies with other neurological diseases (Achard et al., 2006, Lynall et al., 2010, Sanz-Arigita et al., 2010). Additionally, the networks of all individuals had high global Efficiency with broad scale degree distributions.

No significant of HIV status is seen for any of the measures computed (Table 3). However, Modularity, Degree-Variance as well as Assortativity computed from the MCA-GRBF network matrices show trending differences (p < 0.1, uncorrected). On average, the networks of HIV + subject, derived from both methods had decreased Clustering Coefficient (~ 2%), Degree-Variance (~ 12%), and increased Efficiency (~ 3%), though these changes were not statistically significant.

Network hubs based on betweenness centrality

In a preliminary regional analysis, we identified specific regions that serve as “control” regions or “hubs” we computed the mean betweenness centrality for each node across subjects in each group using the MCA-GRBF derived functional network. We defined hubs as specific nodes with Betweeness Centrality values above one standard deviation. Using this criterion both healthy controls and HIV + subjects revealed hub regions, with there being an overlap between them (Fig. 3). Regions of the striatum, specifically bilateral caudate and left thalamus, the left insula and the temporal region appear as hubs in both subject groups. Portions of the parietal lobe appear as hubs only within the control group. However, additional regions across the brain were identified as hubs in the HIV + subjects. When subjects showing symptoms of HAND (ANI and MND group) were analyzed, a further increase in the number of ‘hub’ regions is seen.

Fig. 3

Hub regions as identified based on non-linear connectivity in controls (left), HIV + subjects (center) and subjects showing symptoms of HAND based on NP testing (right) as described in the text. Additional details of the abbreviation labels are provided in the supplementary material. Node sizes in the figure are proportional to the value of Betweenness Centrality. We see that more regions are identified as hubs within the HIV + subjects. There is significant overlap between the hubs of all groups. However, a gradual redistribution of hub characteristics to surrounding regions is observed in the HIV + subjects (center). Furthermore, an additional redistribution of hub regions is seen in subjects showing MND and ANI level of neurodegeneration.

Relationship between imaging metrics and clinical markers

Furthermore, we assessed whether observed changes in brain networks properties were associated with common laboratory plasma markers of HIV infection (CD4 counts and viral load) or the duration of infection. No significant associations were found between such markers and any of the network metrics.

Discussion

This study shows that advanced time-series analysis techniques and graph theory can detect associations between non-linear connectivity derived metrics and cognitive decline occurring during HIV infection. We have tested and shown significant (p < 0.05, FDR corrected) associations between scores of cognitive performance as measured through detailed neuropsychological testing, (the current clinical standard for detecting the presence of HAND) and parameters derived using MCA-GRBF for measuring connectivity (Table 2, Fig. 4). Furthermore, similar analyses based on connectivity measured using a more conventional approach (namely, cross-correlation between low-pass filtered fMRI time-series) are less sensitive in capturing cognitive impairment. Our results suggest that the additional information captured by the MCA-GRBF approach may provide a sensitive measure of CNS injury in subjects with HIV associated neurocognitive disorder.

Fig. 4

Associations between graph-theoretic measures and NP test z-scores. The color code indicates the negative logarithm (to improve visualizations) of p-values after FDR adjustment. Both axes in the figures have been sorted by the overall p-values. It can be seen here that significant association (p < 0.05, − log(p) > 2.99 and ‘*’ in the figure) is observed between multiple characteristics derived from the MCA-GRBF approach (left panel, top left corner), particularly with Executive and Overall z-scores. Some trending associations (p < 0.1, − log(p) > 2.30 and ‘o’ in the figure) as also noted. In contrast, no significant association is seen with correlation-derived metrics (right panel).

The ubiquitous modular structure of the brain appears to be affected during HAND progression as the network property, Modularity was significantly associated with the Executive functioning z-score and also the overall z-score (p < 0.05, FDR corrected). This property captures the prevalence of high-local neighborhood connectivity with reciprocal inter-regional connections. It is possible that the loss or weakening of inter-module connections could manifest as changes in executive functioning. Alterations of inter-regional connectivity reduction have been reported in HIV + subjects across distinct functional networks (Thomas et al., 2013, Wang et al., 2011). It has also been shown, (Jahanshad et al., 2012) in structural studies, that there is a gradual loss of connections within and from the frontal regions that are mainly responsible for executive functioning. Interestingly, changes to the modular structure of the brain have been previously reported in cases of diseases such as Alzheimer's (Chen et al., 2013), schizophrenia (Alexander-Bloch et al., 2014), epilepsy (Vaessen et al., 2011) and similar effects have also been shown to be associated with cognitive decline occurring due to aging (Gallen et al., 2016).

Another network measure, Small-Worldness, was also significantly associated with both the executive functioning as well as the overall z-score. Although not synonymous, Modularity can be considered as closely reflecting the presence of small world architecture (Meunier et al., 2009). Small-Worldness quantifies the ability of the brain to be fully connected yet optimize processing costs thereby providing a high-level description of the brain network. Changes in small world characteristics have been reported in various neurologic diseases (Baggio et al., 2014, Stam et al., 2007). It is interesting to note that the characteristics summarizing two divergent properties of the brain network namely, functional segregation and integration capture information regarding multimodal domains of cognitive processing. We postulate that Modularity and Small-Worldness better represent the ubiquitous integrative nature of the brain network when compared to the other measures used in this study and hence better characterize the complex functions of the brain. We also noted multiple associations with the Executive z-score (Degree Variance, Modularity, Small Worldness as well as Clustering Co-efficient) when analyzing subjects from only the HIV + cohort. These results however were not significant after FDR adjustment, possibly due to reduced variance of cognitive performance and graph measures as well as reduced statistical power in the segregated analysis. Nonetheless we do note that, a general propensity of graph theoretic properties from a non-linear analysis to be strongly associated with Executive functioning which is responsible for integrated cognitive control of behavior. This could also be due to the fact that in this study our analysis was restricted to global measures of brain connectivity computed from a whole brain analysis, hence it is possible that the summary graph theoretic measures were unable capture information pertaining to predominantly unimodal functions (such as motor, learning or memory as measured by NP testing), which could be predominantly localized to specific regions within the brain.

A detailed regional analysis was not performed in this study we however present some preliminary results pertaining to regional changes captured by our method. We analyzed the changes in the hub profile across the subject groups in our study by measuring the Betweenness Centrality (BC), which is defined as the fraction of all shortest paths in the network that pass through a given node. A large body of literature (Achard et al., 2006, Buckner et al., 2009, Crossley et al., 2014, van den Heuvel and Sporns, 2013) has explored this concept, and determined consistent regions of the brain which play a central role in enabling communication across diverse functional and anatomic regions of the brain. In healthy subjects, using MCA-GRBF, we identified the caudate and the thalamus as important hubs of the whole brain network. Within the cortical regions the hubs identified were portions of the frontomedial cortex, superior as well as inferior parietal lobe (Fig. 3). These studies are in line with previous work (Buckner et al., 2009, van den Heuvel and Sporns, 2013) which reported a strong concentration of hubs within the sub-regions of the default mode network. As sub-cortical regions were included in our analysis we notice a predominance of these regions as hubs, which is in line with the knowledge that the basal ganglia are highly interconnected with the cerebral cortex, thalamus and other regions of the brain. It should be noted here that although there exists some overlap, we do expect some divergence in the hubs detected when connectivity is measured using non-linear approaches as in contrast to conventional correlation based metrics. As such regions are critical for heterogeneous functioning as well as efficient co-ordination in the brain network, altered hub characteristics have been linked to the progression of various neurologic and psychiatric disorders (Caeyenberghs et al., 2012, Shu et al., 2016). When the hubs across the different groups were analyzed, we observed that in the HIV + subjects there is a reproducibility of the hubs but a general redistribution of high BC regions from the parietal to the frontal areas of the brain. Furthermore, within the HIV + group, subjects with cognitive impairment (ANI and MND, (Antinori et al., 2007b)) produced a more pronounced redistribution of hub regions. Researchers have proposed the “hub overload and failure” theory (Crossley et al., 2014, Stam, 2014) as a common pathway for progression of various neurological disorders, through which it is hypothesized that, the brain tends to redistribute the network traffic in vulnerable region to other regions in its network hierarchy. This initially tends to rewiring of the network and can eventually lead to hub failure. It is critical to note that changes across different domains were captured via NP testing indicating the possibility that such a route could also be followed during the progression HAND. A similar effect was noted by (Thomas et al., 2015), as they observed a lowering of node centrality values in HIV + individuals within the Parietal and Default mode networks. Our hub analysis results indicate applicability of our method to capture reproducible hubs and study regional differences; we hope to study this aspect in greater detail in future studies.

Comparison between MCA-GRBF and correlation

A key aspect of this study was to investigate whether our method, MCA-GRBF, could provide additional relevant information pertaining to functional connectivity changes, we performed similar analyses using a more traditional cross-correlation to capture connectivity from the low-pass filtered regional time-series. No significant associations were seen between correlation-derived measures and NP testing parameters, suggesting that MCA-GRBF may be more sensitive in capturing changes in network connectivity, than correlation methods. Similar results were obtained in the segregated analyses in each cohort. We believe this specific advantage is related to the fact that MCA-GRBF provides a non-linear approach to measuring functional connectivity while also avoiding the controversial step of low-pass filtering which often leads to disregarding valuable information, induce spurious correlations and even produce discrepancies in results (Davey et al., 2013). It is well established that the brain is a non-linear network, and hence we hypothesize that a non-linear approach, such as MCA-GRBF, would be better suited to capture regional interactions based on BOLD-fMRI.

It is worthwhile to point out here that the for the target disease in our study, the critical question is not necessarily the distinction between HIV + and HIV − individuals but rather identifying the presence of neurodegeneration during HIV infection. HIV + individuals can present with varied scope of neurodegeneration, which is substantially harder for clinicians to diagnose in contrast to detecting HIV infection. Our cohort of HIV infected individuals mostly showed mild cognitive impairment, thus likely reflecting the small changes in functional connectivity. This in combination with a small sample size and appropriate multiple comparisons adjustment reduced the statistical significance of our results at the level of disease (HIV + and HIV −). Interestingly, our results on detecting associations between NP testing derived parameters (the current clinical standard for diagnosing the presence of HAND) and non-linear connectivity-derived measures are promising underlining the potential for development of novel non-invasive imaging biomarkers. However, certain limitations of our study should be mentioned. In our analyses, we have defined regions based on the AAL template, which reduces the brain to 90 regions and may thus constitute an over-simplification of the brain network. We chose the AAL template, because it continues to be widely used in literature, though recent studies have discussed possible advantages of more fine-grained parcellation schemes (Smith et al., 2011). We also acknowledge that the small size of our dataset may not be sufficient to draw conclusions regarding our methods however; our results here do indicate a clear advantage of MCA GRBF over conventional correlation. Future studies with a larger clinical cohort will help assess this advantage further. Another limitation is the cross-sectional nature of our study. We thus cannot draw final conclusions on the applicability of our method as a predictor of HAND progression over time, because this would require a longitudinal study design. It will be important to address both limitations in future research.

In conclusion, in this study we have reported associations between topological network characteristics and cognitive performance scores used for detecting the presence of HIV associated neurocognitive disorder. The non-linear MCA-GRBF method in conjunction with graph-theoretic analysis revealed statistically significant associations with neuropsychological testing results. As a comparison, conventional correlation analysis of low-pass filtered fMRI time-series did not produce significant associations with any clinical parameters. Our results suggest that non-linear analysis of resting-state fMRI can be useful for quantitative evaluation of CNS injury. A longitudinal clinical study needs to be performed in order to further investigate the potential of our method for the development of novel imaging biomarkers for early detection and prediction of disease progression in HAND. The methods proposed here, however can be more generally applicable to the detection of neurodegeneration.