nif:isString
|
-
The proposed method to perform integrated analysis of functional and structural connectivity is illustrated in Fig 1. Each subject requires a DWI, fMRI and structural T1-weighted image. First the T1 image is segmented into a set of regions using FreeSurfer [17], after which structural and functional connectivity maps are extracted. The connectivity maps of each subject are compared uni- and bi-modally for their utility in group-wise analysis as well as identification of connections that contribute to group differences and similarities. The methods are applied to two proof-of-concept studies, one that compares groups of middle-age and elderly subjects, and the second that compares patients with schizophrenia to controls. More details of these steps are provided below.
Figure data removed from full text. Figure identifier and caption: 10.1371/journal.pone.0137484.g001 Schematic of method.Networks are extracted for each individual (upper panel) by extracting brain regions from the T1 image using FreeSurfer, and extracting structural and functional connectivity between these regions from preprocessed DWI and fMRI, respectively. Subsequently, integrated analysis of group differences in structural and functional connectivity is performed (lower panel). Group differences are visualized on a global level using bi-modal distribution plots, on an individual connection level using connectograms and on a mid-level using worm plots, where connections are grouped by location in the brain.
Cortical reconstruction and volumetric segmentation of each subject’s T1 image to produce an ROI atlas was performed using the FreeSurfer image analysis suite [17], available online at http://surfer.nmr.mgh.harvard.edu. The atlas was registered with 12 degrees of freedom to both the DWI’s b0 image and a mean fMRI, also using FSL’s FLIRT algorithm.
The DWI's were preprocessed using AFNI and FSL tools [18]. First they were corrected for motion and Eddy currents and the brain extraction tool was used for skull stripping of the Eddy-corrected DWI. Duplicate b0 volumes were merged into a single volume by taking the median across duplicate images at each voxel. FA and mean diffusivity (MD) images were generated from the skull-stripped DWI using DTIFit. All three images (DWI, FA and MD) were masked using a mask generated from the registered FreeSurfer image, including cortex, white matter and subcortical regions. The network extraction step of the “Statistical Analysis of Minimum cost path based Structural Connectivity” (SAMSCo) framework [4] was used to define structural connectivity because it minimizes the influence of directional uncertainty while finding globally optimal paths. A minimum cost path network (DMC) was extracted, with nodes defined by ROIs in the registered atlas and edges defined as minimum cumulative cost to travel between ROIs. A second metric of structural connectivity, the average fractional anisotropy (DFA) network, was produced by dividing the sum of fractional anisotropy along a minimum cost path by the Euclidean distance of the path. SAMSCo produces two values per connection, one from the seed to target region and the other from the target to seed region. Therefore, the SAMSCo networks were symmetrized by choosing the minimum of the two costs (and corresponding DFA values) to represent the edge between each pair of regions. It is important to note that SAMSCo is a less traditional choice of structural connectivity measure, an aspect which is discussed in more detail in the “Image preprocessing and network generation” subsection of the Discussion section.
Preprocessing of the fMRI's was accomplished using ANFI and FSL tools [18,19]. The fMRI images were slice-timing and motion-corrected, both using AFNI tools. The first four TRs were removed to allow for equilibration and a temporal 0.01Hz high pass filter was applied to remove low frequency signals. Global mean, mean cerebral spinal fluid, and mean white matter time series as well as the six motion parameters resulting from motion correction were regressed out from the image. Finally, an 8mm full-width-half-maximum Gaussian spatial filter [20] was applied to increase the signal to noise ratio. We defined two atlas-based metrics for defining functional connectivity. Typically, functional connectivity studies include a preprocessing step where fMRI images are spatially normalized to a common reference space and blurred to account for differences in structural and functional variability. However, anatomical and functional variation between subjects can be large [20], and thus the deformations involved in spatial normalization may stretch or contract parts of the brain for some subjects more than others. In this study, spatial normalization was not performed. Rather, the FreeSurfer ROIs were registered to the fMRI image for each subject separately, and then used to define connectivity. This reduced the effects of structural variability. The average time series for each ROI was calculated. For the first network type, connectivity between two ROIs was defined as the Pearson correlation coefficient between each ROI's mean time series (fMT). For the other type of connectivity we used partial correlation (fPC), which was calculated for every ROI pair with the influence of all other ROI’s removed according to the formula rij=−pij/piipjj [21], where pij is the (i,j)th element of the inverted covariance matrix from the mean time-series data of all ROIs.
Integrated analysis and visualization of structural and functional networks: An ordinary least squares linear regression model was estimated for each connection including group membership, modeled as a single variable containing zeros and ones, and covariates as variables, along with a constant term to model the mean of the data. The t-statistic and p-values associated with group membership were used to determine whether a connection’s weight differs significantly between groups. The group mean of each connection was also computed for each group of subjects. The t-statistic and group means were visualized using three types of visualizations, illustrated in the lower panel of Fig 1. These include bi-modal network comparison plots of group mean networks and t-statistic networks, worm plots and connectograms. This section provides an explanation of each. All of the visualizations presented in this manuscript can also be used to study continuous variables rather than group differences, which is discussed in more detail in the Discussion section. Bi-modal network comparison plots enable simultaneous evaluation of multiple network types on a global level, both uni- and bi-modally. As seen in Fig 1, they can show distributions of group mean connection weights from two groups simultaneously, or distributions of t-statistics, which illustrates degree of group differences on a global level. Each type of connectivity is represented by a column or row, where rows represent the functional measures and columns represent structural measures. Each connectivity type is evaluated uni-modally in one-dimensional (1D) histograms. The relationship between structure and function is shown in two-dimensional (2D) histograms for each pair of structural/functional network types, which is displayed in the center of the plots. In the case of the group mean networks, colors are used to differentiate between groups, where blue represents one group and red represents the other. In the 2D histograms of group mean networks, magenta indicates group overlap. In the plot of the t-statistics, green is used to represent group differences. Worm plots show the connection t-statistics grouped by the location of the associated ROIs. Even when very few connections are individually identified as significant, the worm plot may show clusters of connections that have significant group differences when considered together (regional shifts in connectivity between multiple ROIs). Significance is thus evaluated on two scales, namely on the level of individual connections and on the level of groups of connections. Worm plots are a derivative of Manhattan plots, which are scatter plots used to show significance of a large number of tests, most of which have very large p-values. They are frequently used in genetic studies [22], where genomic coordinates are on the x-axis and colors of points represent groupings of points. A few key changes to the Manhattan plot were made to make it more suitable to brain connectivity analysis. First of all, connections were clustered by the parts of the brain to which its ROIs belong. Unlike in genetic studies where base-pairs are ordered based on their position within DNA, in the brain there is no intrinsic ordering of connections within each cluster, therefore each cluster was ordered according to t-statistic values. This reordering gives an indication of the distribution of points within each cluster, and enables simultaneous plotting and comparison of multiple network measures. Normally the y-axis of a Manhattan Plot is −log10 (p) where p is the p-value corresponding to a statistical test. We instead used −log10(p)*sign(t)*s0.05,Meff,G, where t is the t-statistic of a given connection, p is the corresponding p-value and sα,Meff,G=log10(α)/log10(αMeff,G) is a scaling factor used to ensure that the reference line corresponding to the adjusted significance threshold, αMeff,G (see below for details of its derivation), is the same for all networks G. sign(t) adds information about whether an observed association is positive or negative. This, combined with ordering connections, resulted in clusters of points with a worm-shape, hence the name Worm Plot. The distance of the mean of each cluster from zero indicates the degree to which the associated parts of the brain play a role in the relationship being tested. Significant connections for each network were also plotted in a connectogram, which was originally used to study brain connectivity by van Horn et al. (2012) [14]. The 2D circular representation in a connectogram provides a cleaner view of brain connectivity than 3D representations. While 3D representations can provide spatial information about ROIs [7,12,13,15,16], viewing them requires projections onto a 2D space. The loss of one dimension results in ROIs appearing spatially closer than they actually are. For example, when viewing a 3D representation axially, the frontal pole and lingual gyrus appear to be located beside each other, when in reality they are on opposite sides of the brain. This also presents difficulties in interpreting region groupings, since groups which are spatially separate will likely appear to overlap in the projected view. The connectogram removes these complications by disregarding 3D spatial information. In the circular representation, ROIs are arranged on the edge of the circle, grouped into the same clusters as in the worm plots. Each cluster is assigned a color, and the concentration of color is adjusted according to node degree. Significant connections are represented as lines between ROIs. The color of each significant connection is determined by the sign of the t-statistic and the opacity indicates the degree of significance. Both worm plots and connectograms require a threshold of test significance. Given that a t-statistic is computed for each connection, some form of multiple testing correction must be used. In the case of a symmetric network G with r regions, the total number of tests is MG = r(r−1)/2. Traditional methods, such as the Bonferroni and Šidák methods, assume independence of tests and calculate a new significance threshold from the total number of tests. However, the connected nature of brain networks implies dependence between connections. As such, traditional methods are overly conservative for connectivity analysis. To address this, we used a method proposed by Li et al (2012) [23], where an “effective number of independent tests” for network G, Meff,G, is calculated. Meff,G replaces MG in the traditional methods. Meff,G was applied to the Bonferroni correction to calculate an adjusted p-value threshold, αMeff,G=α/Meff,G. We used α = 0.05 as the uncorrected threshold.
Two proof-of-concept studies were performed using the analysis described above. The first compared middle-age and elderly subjects from the Rotterdam Scan Study (RSS) [24,25]. The second compared patients with schizophrenia to healthy controls from the Massachusetts General Hospital (MGH) in the MIND Clinical Imaging Consortium data collection [26], which is an open-access multi-site collaborative study of patients with schizophrenia. The original MGH imaging data and subject meta-data is available for download via http://coins.mrn.org/dataexchange. The networks calculated for each subject and the regression results calculated as part of this study are available at http://dx.doi.org/10.5061/dryad.88q04. There were 84 and 82 cortical and subcortical regions chosen for network construction for the MGH and RSS studies, respectively. The cerebellum was excluded in the RSS sample because it was only partially scanned. After registration of the FreeSurfer segmentation to the diffusion image space, which required up-sampling of the data, two regions of one subject in the RSS study did not have any voxels remaining. Regional statistical tests involving these regions were computed by substituting group means for the missing values. See S1 Table for a list of regions used in these studies. In the RSS data, after registration and down-sampling of the FreeSurfer atlas to the fMRI-space, in two regions of one subject no voxels survived the down-sampling (the left entorhinal cortex and the right temporal pole). All analyses did not include the missing values. In the RSS study, non-demented subjects with T1, rs-fMRI and DWI scans were divided into two groups (middle-aged and elderly) based on age. Age ranges were 74 to 95 in the elderly group and 51 to 53 in the middle-aged group. They were selected such that the gender ratio in each group was the same. Subjects with mild cognitive impairment (as determined by the Mini-Mental State Examination [27], MMSE), large infarcts, and gliomas were excluded from the study. The MMSE score is a derived measure of a subject’s global cognitive functioning with a total possible score of 30 [27]. In the MGH study, groups of patients with schizophrenia were compared to controls. Exclusion criteria included an intelligence quotient less than 70, a history of significant head injury and a contraindication for MRI scanning. Only subjects who had DWI, fMRI and pre-computed FreeSurfer segmentations were included. Parental socio-economic status (P-SES) was defined by the Hollingshead-Redlich scale [28], in which a score of one defines the highest status and five defines the lowest. P-SES is often used as an estimate of what the socio-economic status of patients would have been if they did not have schizophrenia. Years of education is defined as the total number of years that a subject pursued elementary, secondary and post-secondary education. In both studies, subjects with more than 4mm of head movement were excluded. All participants gave written informed consent and the original studies were approved by the medical ethics committee of the Erasmus MC, University Medical Center Rotterdam (for the RSS study) and the institutional review board of the Massachusetts General Hospital (for the MGH study). Table 1 describes the demographics of each group for the subjects that were included.
Table data removed from full text. Table identifier and caption: 10.1371/journal.pone.0137484.t001 Subject demographics. * One subject was missing this information, statistics computed over the remaining subjectst Calculated using Welch’s t-testχ Calculated using χ2 test All subjects from the RSS were scanned on the same 1.5 Tesla GE scanner with an 8-channel head coil. A research physician visually inspected each scan to ensure good quality. The T1-weighted protocol included a 3D fast RF spoiled gradient recalled acquisition in steady state with an inversion recovery prepulse sequence (TR = 13.8 ms, TE = 2.8 ms, TI = 400 ms, FOV = 21 × 21 cm2, matrix = 416 × 256, zero-padded to = 512 × 512, flip angle = 20°, NEX = 1, bandwidth = 12.50 kHz) with 96 contiguous slices with slice thickness of 1.6 mm zero-padded to 0.8mm. The final voxel size was 0.49 × 0.49 × 0.8 mm. The DWI protocol utilized a single shot, diffusion-weighted spin echo-planar sequence (TR = 8000 ms, TE = 68.7 ms, FOV = 21 × 21 cm2, matrix = 96 × 64, zero-padded to 256 × 256, NEX = 1) with 36 contiguous slices with 3.5 mm slice thickness. The maximum b-value was 1000 s/mm2 in 25 non-collinear directions. Volumes were also acquired without diffusion weighting (b0). The final voxel size was 0.8 × 0.8 × 3.5 mm. The rs-fMRI protocol included a gradient-echo BOLD sequence, (TR = 2900 ms, TE = 60 ms, matrix = 64 × 64, flip angle = 90°) with 35 contiguous 3.3 mm slices and 160 volumes (total scan length = 464 s). The voxel size was 3.3 × 3.3 × 3.3 mm. The MGH subjects were scanned on a 3 Tesla Siemens Trio scanner with an 8-channel head coil. The imaging parameters of the T1-weighted sequence were TR = 2530 ms, TE = 3.79 ms, TI = 1100 ms, FOV = 16 mm, matrix = 256×256×128 cm, flip angle = 7°, bandwidth = 181, 0.625×0.625 mm voxel size and slice thickness 1.5 mm. The DWI scan parameters were TR = 10,500 ms, TE = 98 ms, NEX = 2, bandwidth = 1342, 64 slices, thickness = 2 mm, 2×2 mm voxel size. The maximum b-value was 1000 s/mm2 in 12 directions. The fMRI parameters were TR = 2000 ms, TE = 30 ms, flip angle = 90°, FOV = 22 cm, bandwidth = 3126Hz/pixel, 27 slices, slice thickness = 4 mm with 1 mm skip, 3.4 mm in plane resolution and 177 volumes (total scan length = 354 s). Three runs of the Sternberg item recognition paradigm (SIRP) [29] were collected. SIRP is a task-based fMRI sequence with a block design that measures working memory. The run with the least amount of motion was selected for connectivity analysis.
Covariates were included in the regression model to control for their influence. In the RSS study covariates included MMSE and gender. In the MGH study covariates included age and gender.
|