Scalable Algorithms for Generating and Analyzing Structural Brain Networks with a
Varying Number of Nodes
Yu Jin1 , Joseph F. JaJa1 , Rong Chen2 , Edward H. Herskovits2
1 Institute

for Advanced Computer Studies

Department of Electrical and Computer Engineering
University of Maryland, College Park, USA
2 Department

of Radiology

arXiv:1609.03893v1 [] 13 Sep 2016

University of Maryland, Baltimore, USA

Abstract
Diffusion Magnetic Resonance Imaging (MRI) exploits the anisotropic diffusion of water molecules in the brain to enable
the estimation of the brain’s anatomical fiber tracts at a relatively high resolution. In particular, tractographic methods
can be used to generate whole-brain anatomical connectivity matrix where each element provides an estimate of the
connectivity strength between the corresponding voxels. Structural brain networks are built using the connectivity
information and a predefined brain parcellation, where the nodes of the network represent the brain regions and the edge
weights capture the connectivity strengths between the corresponding brain regions. This paper introduces a number of
novel scalable methods to generate and analyze structural brain networks with a varying number of nodes. In particular,
we introduce a new parallel algorithm to quickly generate large scale connectivity-based parcellations for which voxels in
a region possess highly similar connectivity patterns to the rest of the regions. We show that the corresponding regional
structural consistency is always superior to randomly generated parcellations over a wide range of parcellation sizes.
Corresponding brain networks with a varying number of nodes are analyzed using standard graph-theorectic measures,
as well as, new measures derived from spectral graph theory. Our results indicate increasingly more statistical power
of brain networks with larger numbers of nodes and the relatively unique shape of the spectral profile of large brain
networks relative to other well-known networks.
Keywords: Connectivity-based brain parcellation, Structural brain network, Spectral clustering, Network analysis,
Spectral analysis

1. Introduction
Non-invasive neuroimaging technologies, combined with
computational techniques, have played a central role in
providing critical insights into the structural and functional organization of the brain. Corresponding brain networks have been extensively studied using different notions
of connectivity such as structural connectivity, based on
synapses between neighboring neurons or fiber tracts between brain regions; functional connectivity based on statistical correlations of activation levels between distinct
regions; or effective connectivity that captures causal interactions that can be inferred from network perturbations
or time-series analyses. These networks have been studied using segregation measures that refer to the degree
to which a network’s nodes form separate cliques or clusters, and integration measures that refer to the capacity
of the network as a whole to become interconnected and
exchange information [1].
This work focuses on structural brain networks derived
Preprint submitted to NeuroImage

from Diffusion Tensor Imaging (DTI) data. DTI allows one
to observe the diffusion process of water molecules in brain
tissues under magnetic fields with different strengths and
along different directions [2, 3, 4]. Modern tractographic
methods can then be used to estimate the anatomical connectivity information. Using this connectivity information and a set of predefined regions of interests (ROIs) or
whole brain parcellations, structural brain networks can be
built and analyzed. A significant number of efforts have
been devoted to study the topological properties of structural brain networks using graph-theoretic measures such
as characteristic path length, clustering coefficient, local
and global efficiency, degree distribution, betweenness centrality etc [5, 2, 6]. These studies have resulted in a number of interesting findings regarding brain networks such
as scale-freeness, small-world property, and modular organization. Moreover, researchers have explored the topological differences of the networks corresponding to distinct
population groups [7, 8, 9]. It was observed that subjects
with certain brain disorders seem to have structural abSeptember 14, 2016

• Development of a fast parallel algorithm to generate
whole-brain connectivity-based parcellations with a
very large number of regions. The resulting regions
are shown to have a very good degree of structural
consistency in terms of the corresponding brain networks.

normalities in terms of local and global graph properties.
Given the early promising results, a number of tutorials
and tools have been developed to facilitate brain network
analysis [5, 10].
As described in the literature, structural brain networks have been constructed based on predefined parcellations that typically use traditional anatomical brain atlases, such as the Automated Anatomical Labeling (AAL)
atlas [11, 12] and the Jülich histological atlas [13], or some
random parcellation[6]. Researchers have also developed
parcellation algorithms that are based on well-known clustering methods such as Gaussian Mixture Model (GMM),
k-means clustering, hierarchical clustering [14]; these algorithms partition the brain mask into many spatially contiguous regions. These studies have primarily focused on
specific anatomical regions of the brain such as the human inferior parietal cortex complex (IPCC) [15], the lateral parietal cortex [16], the temporoparietal junction area
(TPJ) [17], the dorsal frontal cortex [18], the ventral frontal
cortex [19], cingulate and orbitofrontal cortex [20], and
Broca’s areas [21]. Most of these methods are not easy
to generalize to determine whole-brain brain parcellations
while enforcing the constraint of connectivity-based homogeneous regions. Note that the choice of the parcellation
plays a critical role in determining the quality of the corresponding brain networks. Traditional atlases and random
parcellations do not make use of the connectivity information that can be derived from tractographic methods,
which may alter the statistical effectiveness of subsequent
structural network analysis. Another important parameter that received relatively little attention in the literature
is the number of nodes in the brain network. In particular, it would be interesting to figure out how the number
of nodes affects various network properties and whether
larger brain networks capture additional structural patterns than smaller ones.
The graph theoretical analysis is usually performed on
sparse brain networks derived from the original relatively
dense brain networks [2]. The original brain networks may
contain weak connections that are eliminated by a thresholding process typically based a single threshold value. But
this threshold value affects the structure and sparsity of
the brain networks and hence may heavily impact the corresponding graph properties. In fact, often the threshold
value is selected to enforce certain sparsity of the graph
without fully understanding the impact on the corresponding topological properties of the derived netowrks. During the thresholding (and binarizing) process a significant
amount of connectivity information is lost. This paper
explores the use spectral graph theory to be introduced
later for analyzing brain networks. Spectral analysis is
performed on the original weighted brain network and can
capture almost all the topological invariants of the network. No threstholding process is needed, and the corresponding spectral analysis seems to be much more robust
as we will see later.
The main contributions of this paper are the following.

• Analysis of structural brain networks for distinct population groups using a varying number of nodes, up
to 500 nodes, which indicates richer structural patterns with larger numbers of nodes.
• Preliminary indication of the robustness and the power
of network parameters derived form spectral graph
theory when compared to the traditional graph-theoretic
measures. In particular, the shape of the spectrum
of such networks seems to be unique among the spectrums of other well-known networks.
The rest of the paper is organized as follows. Section 2
introduces the dataset used in this study as well as the software tools used to generate the connectivity information.
Section 3 briefly describes the GPU-CPU heterogeneous
system and necessary software libraries used for our parallel implementation. Section 4 provides a brief overview
of the parallel parcellation algorithm while Section 5 covers the evaluation of the connectivity-based parcellations
and the network analysis using both the traditional graphtheoretic measures and new parameters derived from spectral graph theory. We conclude the work in Section 6.
2. Materials
2.1. Diffusion MRI Dataset
The diffusion MRI dataset used in this research is taken
from the publicly available Nathan Kline Institute (NKI)/
Rockland dataset1 . The neuroimage sample consists of
data for individuals whose ages range from 4 to 85. The
demographics of our sample are shown in Table 1. A
comprehensive psychiatric interviews and behavioral assessments are included to help explore the brain-behavior
relationship. The diffusion MRI was performed using a
SIEMENS MAGNETOM TrioTim syngo MR B15 system.
The high-angular resolution diffusion imaging protocol was
used to assess white matter integrity as measured by fractional anisotropy. Diffusion tensor data were collected
using a single- shot, echo-planar, single refocusing spinecho, T2-weighted sequence with a spatial resolution of
2.0×2.0×2.0mm. The sequence parameters were:
TE/TR=91/10000ms, FOV=256mm, axial slice orientation with 58 slices, 64 isotropically distributed diffusion
weighted directions, two diffusion weighting values (b=0
and 1000s/mm2 ) and six b=0 images. These parameters
were calculated using an optimization technique that maximizes the contrast to noise ratio for FA measurements.
1 http://fcon

2

1000.projects.nitrc.org/indi/pro/nki.html

Table 1: Demographics of NKI Sample

Age Group
Female
Male
Total

4-9
4
5
9

10-19
18
16
34

20-29
19
24
43

30-39
8
9
17

40-49
5
26
31

50-59
6
8
14

60-69
9
5
14

70-85
10
4
14

Total
79
97
176

For each subject, the image data consists of 76 volumes of
3D images of dimensions 128×128×53, each voxel representing 2.0mm×2.0mm×2.0mm brain volume.
3. Computational Environment
3.1. CPU-GPU Heterogeneous Platform
We use a heterogeneous CPU-GPU platform to generate and analyze our structural brain networks. Standard
heterogeneous platforms consist of a multi-core CPU and
one or more general-purpose GPU (GPGPU). The CPU
and GPU communicate through the Peripheral Component Interconnect Express (PCIe) bus. The multi-core
CPU is programmed as a shared memory system with a
large main memory and multiple levels of cache. Most
existing software packages are implemented based on single or multi-core CPU architecture. On the other hand,
GPUs typically have thousands of cores that can execute
data parallel operations extremely fast. They have been
used extensively for all types of computations especially
scientific computations. However programming a GPU is
typically a non-trivial task. Since the CPUs and GPUs
have somewhat complementary strengths, modern applications try to combine them in single platforms [22]. Hybrid
computing environments, which collaboratively combine
the computational advantages of GPUs and CPUs, further boost the overall performance [23, 24, 25, 26, 27, 28].
However such a performance gain is typically achieved at
the expense of a significant additional programming complexity. The details of the platform used in this study are
shown in Table 2.

Figure 1: CUDA Program Model

implement our brain parcellation scheme on CPU-GPU
heterogeneous platforms.
3.2.1. CUDA Library
The Compute Unified Device Architecture (CUDA) is
a general-purpose multithreaded programming model that
utilizes the large number of GPU cores to solve complex
computational problems. The CUDA programming model
assumes a heterogeneous system with a host CPU and
one or more GPUs as co-processors. Each GPU has an
array of Streaming Multiprocessors (SM ), each of which
has a number of Streaming Processors (SP ) that execute
instructions concurrently. The parallel computation on
GPU is invoked by calling customized kernel functions using thousands of threads. The kernel function is executed
by blocks of threads independently. Each block of threads
can be scheduled on any Streaming Multiprocessor (SP)
as shown in Figure 1. The kernel function takes as parameters the number of blocks and the number of threads
within a block.
In addition, we use the efficient BLAS libraries for both
sparse2 and dense3 matrix computations. Our implementation also relies on the Thrust library, which resembles
the C++ Standard Template Library (STL) that significantly improves productivity.

Table 2: System Details

CPU Model
CPU Cores
DRAM Size
GPU Model
Device Memory Size
SMs and SPs
Compute Capability
CUDA SDK
PCIe Bus

Intel Xeon E5-2690
8
128GB
Tesla K20c
5GB GDDR5
13 and 192
3.5
7.5
PCIe x16 Gen2

3.2. Library Dependencies
2 http://docs.nvidia.com/cuda/cusparse/

Our brain parcellation scheme is based on a fast parallel implementation of the spectral clustering algorithm
described in [29]. We leverage the following libraries to

3 http://docs.nvidia.com/cuda/cublas/

3

Figure 2: Overview of the Generation and Analysis of Structural Brain Networks

3.3. ARPACK Software
ARPACK is a software package designed to solve largescale eigenvalue problems [30]. It contains highly optimized Fortran subroutines that solve symmetric, non-symmetric
and generalized eigenvalue problems which are widely used
in modern scientific software packages such as Matlab and
Python scientific packages. ARPACK is based on the
Implicitly Restarted Arnoldi Method (IRAM) with nonFigure 4: Illustration of connectivity profile for each voxel
trivial numerical optimization techniques [31, 32]. The
eigenvalue problem is efficiently solved by collaboratively
target regions. The JHU DTI-based white matter atlas4 is
combining the interfaces of ARPACK and cuSPARSE lispecified as the seed region. The AAL mask is specified as
brary.
the target region, which is also considered as our region of
interests to be parcellated. We generate 500 streamlines
from every voxel in a seed region. These streamlines are
4. Methods
propagated following the cross fiber distribution derived
An overview of the workflow for generating and anafrom the voxel-level diffusion tensors. Curvature threshold
lyzing structural brain networks is illustrated in Figure 2.
is enforced to eliminate unqualified streamlines. The dis. The DTI preprocessing step extracts the connectivity
tance correction option is set to correct for the fact that
information represented as a large-scale connectivity mathe distribution drops as travel distance increases. The
trix by using probabilistic tractography. A scalable paraltractography output is a structural connectivity network
lel algorithm is then used to generate parcellatioins with
modeled as a large weighted graph where each node is a
up to 500 regions [29, 33]. Structural brain networks are
voxel in the target region space and each edge weight corconstructed from the parcellations and the connectivity
responds to relative connectivity strength in terms of the
information. Finally we apply both traditional graphnumber of streamlines connecting the corresponding pair
theoretical methods and new methods based on spectral
of voxels. Figure 3 shows the seed and target regions as
grpah theory to analyze the brain networks. The details
well as the tractography result of a subject.
of each step are described below.
4.2. Brain Parcellation Scheme
4.1. Dataset Preprocessing
The brain parcellation scheme takes as input the probabilistic tractography results represented as a large con4.1.1. Nonlinear Registration
nectivity matrix. The number of voxels in the AAL mask
The diffusion images of all subjects are registered to
is 155, 794 and the connectivity matrix is a sparse matrix
Montreal Neurological Institute (MNI) standard space usof size 155, 794 × 155, 794. Given a positive integer value
ing the nonlinear registration package FNIRT in FSL softk, the brain parcellation problem is to segment the brain’s
ware [34]. The nonlinear registration process generates the
grey matter into k spatially contiguous regions, such that
warping coefficients that balance the similarity between
the connectivity profiles of the voxels in each region are
the diffusion image and the standard MNI152 image, and
as similar as possible. Moreover, we expect the parcellathe smoothness of the warping coefficients.
tions to be stable and reproducible, as well as, consistent
among members of a structurally homogeneous population
4.1.2. Probabilistic Tractography
sample..
The probabilistic tractography method is one of the
We next introduce our notion of a connectivity profile
modern tractography methods which has been used extenfollowed
by an overview of the methods used in our parsively to extract anatomical connectivity information. For
cellation
algorithm.
each voxel, the diffusion tensor is built to model the fiber
distribution through the BEDPOSTX package in FSL [35].
The probabilistic tractography is then processed through
the probabilistic tractography toolbox in FSL [36] by sending out streamlines from seed regions and propagate through

4.2.1. Connectivity profile
For each voxel, the connectivity profile is the signature
that discriminates a voxel from the rest of voxels based on
4 http://cmrm.med.jhmi.edu/

4

Figure 3: Left: Seed region which is JHU white matter atlas. Middle: target region which is the AAL mask. Right: the result of probabilistic
tractography which models the distribution of neuron fiber bundles where the intensity of each voxel represents the number of streamlines
passing through that voxel. All figures are imaged in the axial plane

connectivity information. Parcellations are built by clustering voxels with similar connectivity profiles together.
As shown in Figure 4, the connectivity profile of a voxel
is defined as an array of weights, where each weight represents the connectivity strength from that voxel to the
voxels within each region defined by a given brain segmentation. As shown in [33], our parcellation algorithm
generates almost identical parcellations regardless of the
initial brain segmentation used to define the connectivity
profiles as long as the segmentation has a sufficient number
of regions.
Our brain parcellation scheme is an iterative algorithm
that is briefly summarized in Algorithm 1, and is explained
in detail in [33]. Each major step of the brain parcellation
scheme is described next.

Figure 5: At most 32 neighbors of a voxel are within a sphere of
radius r = 2.

The standard spectral clustering method can be described as follows, where W ∈ Rn×n is the weight matrix
associated with the spatial similarity graph and k is the
number of desired regions.

4.2.2. Spatial similarity graph
The spatial-constraint similarity graph is formed using
spatial adjacency and the connectivity profiles as follows.
The voxels define the nodes of the graph. Two nodes are
connected by an edge if and only if the corresponding voxels are within a spatial distance of radius r. In this study,
we use r = 2 and hence the number of neighbors of each
node is at most 32 as shown in Figure 5. Each edge is
weighted by the correlation coefficient between the connectivity profiles of its end points.

• Compute the normalized Laplacian matrix L = D −
W , where D is P
the diagonal matrix such that each
n
element Di,i = j=1 Wi,j .
• Compute the k eigenvectors of the normalized Laplacian D−1 L (= I−D−1 W ) corresponding to the smallest k eigenvalues.
• Apply the k-means clustering algorithm on the rows
of the eigenvectors to obtain the final clusters.

4.2.3. Spectral clustering
The spectral clustering algorithm is initially proposed
to solve minimum graph-cut problems. Specifically, the
algorithm partitions a similarity graph into several subgraphs with the objective of minimizing the total weight
of the edges connecting the subgraphs subject to a constraint on the relative sizes of the subgraphs [37, 38, 39].
In our case, the subgraphs induce a spatially contiguous
segmentation of the AAL mask. The algorithm results in a
solution where the voxels within the same region have similar connectivity profiles and voxels across different regions
have dissimilar connectivity profiles.

The spectral clustering algorithm is very computationally expensive especially when we are dealing with a large
number of voxels n and a large number k of regions. Therefore we develop high performance computing techniques to
speed up the algorithm as described in the next section.
4.3. Parallel spectral clustering
A detailed description of our parallel spectral clustering
algorithm and its application to a wide variety of datasets
has been presented in [29]. Here we provide a high level
5

Algorithm 1 Connectivity-based Brain Parcellation Scheme
1. Generate the probabilistic tractography results from diffusion MRI data.
2. Initialize a parcellation as a random spatial brain segmentation, to be used to define the connectivity profile of each voxel.
Repeat
3. Define the connectivity profiles using the current brain parcellation.
4. update the spatial similarity graph using the connectivity profiles.
5. Apply spectral clustering algorithm on the current spatial similarity graph to generate the brain parcellation.
6. Measure the similarity between the new parcellation and the previous parcellation used to define connectivity profiles.
Until the similarity measurement exceeds some threshold.
7. Return the parcellation result.

4.3.2. Parallel eigensolver
The most computational expensive part of the spectral
clustering algorithm is the determination of the k eigenvectors of a large real-valued symmetric matrix corresponding
to the smallest k eigenvalues.
We combine an important feature of the ARPACK
software named reverse communication interfaces and the
CUDA-based sparse matrix-vector multiplication. The reverse communication interfaces are a set of CPU-based
interfaces that facilitate the process of solving large-scale
eigenvalue problems. For each iteration, one of the interfaces TakeStep() provides a n-length vector used as the
input, and the output of sparse matrix-vector multiplication is provided back to the interface Converge() to determine the convergence. In our implementation, the matrixvector multiplication is performed on the GPU side using
cuSPARSE library. For each iteration, the input vector is
transferred from the CPU to the GPU and the output vector is transfered back to the interface. The procedures of
parallel eigensolver are described in Figure 6. A detailed
description can be found in [29].

overview of the implementation on a heterogeneous CPUGPU platform.
4.3.1. Parallel graph construction
We start by computing a sparse representation of the
similarity graph on the CPU side using edge lists. Such
representation is only computed once because all subjects
are registered from the individual local space to the standard space, and hence will have the same patterns of spatial affinity. Each element < i, j > on the edge list is
a single edge that connects two nodes indexed by i and
j within a spatial distance r ≤ 2 . As the connectivity
profiles are known, the computation of edge weights is a
data-parallel task which is computed on the GPU with
each thread executing on a small set of the edges.

4.3.3. Parallel k-means clustering
Modern scientific software packages and toolboxes such
as Matlab5 and scikit-learn provide serial and parallel versions of the k-means clustering algorithm. However their
CPU-based implementations are inefficient for large-scale
problems, especially when the number of clusters k is very
large. We briefly present a highly efficient implementation
based on an open-source project6 which efficiently utilizes
the Thrust and CUDA libraries and achieve very significant speedups.
The eigenvectors corresponding to the smallest k eigenvalues are represented as V ∈ Rn×k , which initially reside
on the CPU memory, where n is the number of voxels and
k is the desired number of regions. The procedure is summarized in Algorithm 2.
4.4. Construction of Structural Brain Networks
The brain networks in our study are built based on the
connectivity-based brain parcellations as well as the con-

Figure 6: A fast implementation to solve the smallest k eigenvalue
problem

5 http://www.mathworks.com/help/stats/kmeans.html
6 https://github.com/bryancatanzaro/kmeans

6

sonable measure for the relative connectivity strength between regions.
Our second definition normalizes the edge weights considering the strength of within-region connections. We first
define the strength of the inter-region and within-region
connections of Ri and Rj as follows.
X
W (va , vb )
(2)
Ci,i =

Algorithm 2 Parallel K-means Clustering Algorithm
1. Transfer the eigenvectors V ∈ Rn×d from the CPU
to the GPU.
2. Randomly select k rows of V as the centroids of
the k clusters stored in C ∈ Rk×d
3. While (the centroids change) do
Compute the pairwise distances S ∈ Rn×k between
rows of V and the centroids.
Update the new label of each node.
Compute the new centroids of the clusters.
4. Transfer the labeling result from GPU to CPU.

va ∈Ri ,vb ∈Ri

X

Cj,j =

W (va , vb )

(3)

W (va , vb )

(4)

va ∈Rj ,vb ∈Rj

X

Ci,j =

va ∈Ri ,vb ∈Rj

nectivity information revealed by probabilistic tractography. The network nodes correspond to the regions in the
parcellation. We build and analyze structural brain networks at different scales which correspond to parcellations
with different numbers k of regions.
Given a parcellation whose regions define the nodes of
the brain network, there is no consensus in the literature
on how the edge weights should be defined. At a high
level, the edge weight reflects the number of streamlines
connecting the corresponding endpoints [6, 40, 41, 9], or
the probability that the two endpoints can be reached from
each other [8, 5]. In some cases, the weights are normalized considering the region volumes as well as the lengths
of streamline trajectories. The key problem with defining
edge weights is that we need to accurately estimate the
region-level connectivity strengths from the probabilistic
tractography results that only involve voxel-level connectivity information. Note that the voxel-level connectivity
information corresponds to an estimate of the number of
streamlines connecting the two voxels, independently of
the intermediate voxels. Many researchers such as [8, 4]
compute the region-level connections by accumulating the
number of streamlines between all voxels within the seed
region and all voxels within the target region. Since some
streamlines may pass through many voxels within the same
region, such estimate will likely lead to a significant overestimation. Based on this observation, we introduce two
notions of region-level connectivity strength W (Ri , Rj ) between two region Ri and Rj . The first notion is defined as
follows.
W (Ri , Rj ) =

max

va ∈Ri ,vb ∈Rj

W (va , vb )

The normalized edge weights are defined by
Ci,j
W (Ri , Rj ) = p
Ci,i Cj,j

(5)

By the Cauchy-Schwarz inequality, the normalized edge
weights are always less or equal to 1. When most paths
through region Ri cross region Rj , the edge weight is close
to 1. Otherwise when the intra-connections within regions
are strong but the inter-connections between regions are
weak, the probability is close to 0.
Since both definitions yield similar results in our subsequent analysis, we use the lower bound of the streamlines
to define the edge weights of our structural brain networks
for the rest of the paper.
4.5. Network Analysis
4.5.1. Graph Theoretical Analysis
As described in section 4.4, the nodes of a structural
brain network correspond to the brain regions defined by
a parcellation, and the edge weights can be defined either
by Equation (1) or (5). Since the results are similar with
either definition of edge weights, we will adopt Equation
(1) for the rest of the paper. In the next few sections,
we present two types of network analysis. The first deals
with the classical graph theoretic measures, and is coupled with an exploration of the relationship between the
graph theoretical measures and the scale of brain networks.
The study is conducted relative to two pairs of population
groups that are expected to possess significant structural
differences. The second type of analysis is based on spectral graph theory and a number of parameters derived from
the spectra of structural brain networks. We evaluate the
robustness of the two types of network analysis in addition to exploring their statistical power relative to distinct
population groups.
Many global and local graph-theoretic measures have
been proposed to study different types of networks including social, biological, and social media networks. Some of
the studies that are related to our work appear in [5, 10].
Here we focus on a relatively few of the most frequently
used graph theoretic measures to study brain networks.

(1)

where W (va , vb ) represents the number of streamlines connecting the two voxels as generated by the tractographic
results.
It is easy to see that W (Ri , Rj ) is a lower bound of the
total number of streamlines connecting the corresponding
regions. For parcellations with a large number of regions,
we expect W (Ri , Rj ) to be reasonably close to the true
value since the number of voxels within each region will
be relatively small. Overall we believe that this is a rea7

• For each node, we want to eliminate edges with small
weights. Since the total number of connections is significantly different across subjects, we first normalize
the associated edge weights such that the sum of outgoing edge weights is equal to 1 and then eliminate
the edges whose value are below a certain threshold
. In our study, we choose  = 0.01.

We note that the graph measures we are about to introduce are for unweighted and undirected graphs U where
ui,j = 1 when there is an edge between nodes i and j, and
ui,j = 0 otherwise. The number of nodes is N .
• Characteristic path length (CPL): The characteristic
path length tries to capture the network integration
and is computed as the average of the shortest paths
between all pairs of vertices.
CP L =

1
N (N − 1)

X

d(i, j)

Wi,j
wi,j = P
j Wi,j

(6)
4.5.2. Spectral Analysis
We introduce some basic concepts from spectral graph
theory, which we believe offers a more robust foundation to
analyze and characterize brain networks. In essence, spectral graph theory explores the fundamental relationship
between graphs and the spectrum of matrices associated
with these graphs, typically either the adjacency matrix,
or the Laplacian matrix. There are many tutorials and
several textbooks that have been written about spectral
graph theory. See for example [43].
In our case, we focus on the normalized Laplacian matrix whose eigenvalues are closely related to almost all major topological invariants of the graph. The normalized
Laplacian of a graph (or network) was introduced earlier
in this paper when we gave an overview of the spectral
clustering algorithm. It can be defined using the following
equation:

i,j,i6=j

where d(i, j) is the shortest path between node i and
j.
• Global efficiency (Eglobal ). The global efficiency of
a graph is the average of the inverse of the shortest
paths between all pairs of vertices and hence tries to
capture how well pairs of nodes are connected. [42].
Eglobal =

X
1
1
N (N − 1)
d(i, j)

(7)

i,j,i6=j

• Clustering coefficient. The clustering coefficient tries
to capture graph separation and is defined as the average of the local clustering coefficient at each node.
C=

1 X
Ci
N

(8)
I − D−1 W

where the local clustering coefficient at node i is defined as
2Γi
Ci =
(9)
degi (degi − 1)

1X
ui,j ui,h uj,h
2

(13)

where W is the weight matrix of the graph and D is the
diagonal matrix
such that each diagonal entry is defined
P
by Di,i = j Wi,j .
The (normalized) spectrum of the graph is defined by
the sequence of the eigenvalues of the normalized graph
Laplacian ordered as follows: λ1 = 0 ≤ λ2 ≤ ... ≤ λn . It
can be shown that all eigenvalues are between 0 and 2 and
that the largest eigenvalue λn is equal to 2 if and only if
the graph is bipartite [43]. Moreover, the second smallest
eigenvalue λ2 is of special importance and has been studied extensively in the literature. For example, it turns out
that λ2 is related to the convergence rate of random walks
on the graph, ”chemical conductance,” Cheeger’s constant
(roughly the size of the sparset cut in the graph), graph
diameter, and graph quasi-randomness [43]. Other eigenvalues also play an important role in capturing the graph
structural properties.
Compared with traditional graph theoretic measures,
spectral analysis seems to offer the following advantages,

Γi is the number of triangles around node i,
Γi =

(12)

(10)

j,h

Note that degi > 0 in our case since isolated nodes
are removed from the network.
• Sparsity ratio. The sparsity ratio is defined as the
ratio of the number of edges over the number of possible edges between nodes, that is,
P
i,j ui,j
Sparsity =
(11)
k(k − 1)
A more detailed description of graph theoretic measures commonly used to analyze brain networks can be
found in [10]. We preprocess our weighted structural brain
networks into sparse unweighted graphs as follows [5],

• The eigenvalues of the graph Laplacian can be computed directly on the weighted structural brain network, which can be shown to be robust to small
noise.

• Remove self-loops whenever they exist.

• The entire spectrum seems to capture almost all the
invariant structural properties of the network.
8

Figure 7: Example of connectivity-based brain parcellations at multiple scales. (Random coloring; axial view)

We will also illustrate the robustness of the graph spectrum relative to graph theoretic measures.

process from local space to standard space; these voxels
are eliminated before evaluating regional consistency.

5. Results

k
1X
1
k i=1 |Ri |(|Ri | − 1)

5.1. Parcellation Evaluation

X

Corr(Pva , Pvb )

(14)

va ,vb ∈Ci ,vb 6=va

We start by shedding some light on the quality of our
parcellations as the number of regions grows. The connectivity- The regional consistency measures the homogeneity in
terms of the correlation coefficients of the connectivity probased brain parcellations with different number of regions
files averaged over all regions. Figures 8 show the average
k are shown in Figure 7. We compare our connectivityregional consistency of the connectivity-based brain parbased parcellations against randomly generated parcellacellations and random parcellations at multiple scales as
tions which are commonly used in many connectome studwell as AAL atlas over five subjects where the connectivies. The randomly generated parcellations are a set of
ity profiles are defined based on the corresponding parrandom spatial segmentations of AAL mask which are gencellations. The initial segmentation used to define the
erated in a similar way as in [6].
connectivity profiles is the AAL atlas. For all subjects,
We evaluate the regional consistency of any parcellathe connectivity-based parcellations consistently possess
tions as the average of the correlation coefficient of connechigher structural regional consistency than the random
tivity profiles between all pairs of voxels within the same
parcellations (p < 10−10 ) as well as AAL atlas, which
region, which are defined in equation (14). Note that Pva
suggests that connectivity-based brain parcellations have
is the connectivity profile of voxel va defined by the parcelbetter voxel-wise homogeneity that are more suitable for
lation and |Ri | is the number of voxels within the region
constructing structural brain networks. Moreover, the avRi . We justify this notion of connectivity profiles as folerage regional consistency slightly decreases as the number
lows. The goal of our pacellation is to build structural
of regions increases mainly because the connectivity probrain networks, and hence the voxels captured by a node
files of fine-grained parcellations capture more connectivity
should have very similar connectivity patterns to the rest
information than the coarse-grained connectivity profiles.
of the nodes of the network. Otherwise we won’t be able
to compress all the voxels in a region into a single node.
We also note that some voxels have no or very few connections to other regions due to the nonlinear registration
9

Figure 8: Regional consistency of connectivity-based brain parcellations and random parcellations.

Table 3: Group Information

Group Number
Age Range
Female
Male
Total

Comparison 1
Group 1 Group 2
4-10
60-80
7
17
5
7
12
24

decreases only very slightly with increasing k. These trends
are consistent with prior results in the literature such as
[6]. We will shed more light on the effect of sparsity ratio
on various graph-theoretic measures later in this paper.
We now examine the statistical differences of our graph
theoretic measures between two sets of pairs of population
groups that are expected to possess different structural
patterns. The two sets are defined in Table 3. Statistical
differences are measured using the p-values of the twosample t-test, over a range of network sizes.
Table 5 includes the p-values for each group pair and
when the number of nodes varies from k = 50 up to k =
500.
We note that the characteristic path length and global
efficiency are statistically different for each pair of population groups, regardless of the number of nodes in the network. However the clustering coefficient becomes statistically discriminative for each pair of the population groups
only when the number of regions becomes relatively large,
at least k ≥ 200. This illustrates the fact that larger brain
networks can capture more detailed structural differences
than smaller networks.

Table 4: Graph Measurements of Structural Brain Networks

k
50
60
70
80
90
100
200
300
400
500

CPL
1.8641
1.9610
2.0375
2.0929
2.1102
2.2181
2.6503
2.9978
3.2075
3.4294

Eglobal
0.6159
0.5853
0.5652
0.5488
0.5394
0.5182
0.4326
0.3829
0.3562
0.3328

C
0.5253
0.4917
0.4711
0.4581
0.4527
0.4512
0.4236
0.4268
0.4214
0.4209

Comparison 2
Group 1 Group 2
11-20
60-80
17
17
19
7
36
24

Sparsity
0.2796
0.2345
0.2118
0.1892
0.1708
0.1560
0.0819
0.0541
0.0411
0.0333

5.2. Network Analysis
5.2.1. Graph Theoretical Analysis
Table 4 shows the relationship between the graph theoretic measures introduced above and the number of regions. The global efficiency decreases significantly with
the number of regions while the characteristic path length
increases significantly. In addition, the sparsity ratio decreases quite significantly while the clustering coefficient
10

Table 5: Discriminative Power of Graph Measurements (in terms of p-value of the two-sample t-test)

k
50
60
70
80
90
100
200
300
400
500

CPL
4.0466 × 10−6
4.0429 × 10−7
5.5281 × 10−6
9.5953 × 10−6
6.5179 × 10−6
6.0881 × 10−6
3.0427 × 10−7
5.3190 × 10−5
4.7068 × 10−5
1.4848 × 10−5

Comparison 1
Eglobal
5.2804 × 10−6
1.3503 × 10−6
5.8380 × 10−6
1.1765 × 10−5
1.0460 × 10−5
8.7734 × 10−6
4.8142 × 10−7
7.0366 × 10−5
3.0236 × 10−5
2.1470 × 10−5

C
0.2965
0.8713
0.8108
0.3231
0.0934
0.0351
7.7821 × 10−7
3.2108 × 10−8
1.2254 × 10−6
9.5681 × 10−6

CPL
0.0078
7.8283 × 10−4
5.6084 × 10−4
0.0020
5.1846 × 10−4
3.3375 × 10−4
2.2048 × 10−4
6.5096 × 10−4
7.0862 × 10−4
0.0015

Comparison 2
Eglobal
0.0038
5.6824 × 10−4
6.7248 × 10−4
0.0031
7.2236 × 10−4
5.1019 × 10−4
8.7778 × 10−4
0.0042
0.0097
0.0204

C
0.8352
0.6866
0.2981
0.4535
0.5989
0.2284
3.4852 × 10−4
4.7756 × 10−5
1.5504 × 10−6
6.1970 × 10−5

Table 6: Discriminative Power of Spectral Measurements (in terms of p-value from the two-sample t-test)

k
50
60
70
80
90
100
200
300
400
500

Comparison 1
Modularity
λ2
0.6033
0.0012
0.0223
1.9278 × 10−4
0.0049
1.1783 × 10−4
0.0014
1.0335 × 10−4
0.0010
1.6525 × 10−4
0.0024
1.1762 × 10−4
−5
7.6224 × 10
1.0702 × 10−6
1.0352 × 10−4 1.5826 × 10−6
4.6564 × 10−4 3.1102 × 10−6
3.1176 × 10−4 4.5524 × 10−6

Comparison 2
Modularity
λ2
0.0146
0.0094
0.0987
5.3686 × 10−5
0.0282
2.3541 × 10−4
0.0117
8.0141 × 10−5
0.0238
3.7983 × 10−4
0.0026
5.7496 × 10−5
0.0281
9.0932 × 10−6
0.0275
3.3413 × 10−4
0.0439
3.2788 × 10−4
0.1043
0.0013

Table 7: Prediction Accuracy with Graph Topological Measurements and Spectral Properties of Structural Brain Networks

Feature selection
CPL+Eglobal + C
λ2
λ2 + λ3

k = 500
Accuracy (Group 1) Accuracy (Group 2)
0.8056
0.7167
0.8056
0.7333
0.8333
0.7667

Since the graph-theoretic measures represent discrete
aggregated network measurements, they lack the ability to
capture network characteristics in detail, which motivates
the analysis based on the entire spectrum of the network
of each subject as described next.

k = 90
Accuracy (Group 1) Accuracy (Group 2)
0.7778
0.6500
0.7778
0.7000
0.7778
0.7333

larger networks are further and further away from being
bipartite.
We note that spectral analysis has been recently used
to analyze different types of networks such as social networks and biological networks [44]. Most related to our
work are the unweighted neural networks of macaque, cat
and Caenorhabditis elegans [45].
We introduce two global parameters based on the spectrum of a network.

5.2.2. Spectral Analysis
Figures 9 -12 show the spectral plots (histograms of the
normalized Laplacian eigenvalues) of the weighted structural brain networks of two subjects, each with two distinct
values of the number of nodes. Note that the spectral plots
share very similar shapes, which hold over the subjects in
our data set especially for larger networks. Note also that
the second eigenvalue decreases with the size of the network and larger fractions of smaller eigenvalues appear
with larger networks (which indicate more modular structures). Note also that the spectrum of larger networks is
further and further away from 2, which implies that the

• The second smallest eigenvalue λ2 . As mentioned
before, this parameter seems to capture a number of
critical topological properties of the network.
• Modularity. The number of small eigenvalues in the
spectral plot reflects the extent of structural modularity in the network. In this paper, we define the
modularity as the number of eigenvalues less than a
certain threshold, say γ = 0.3.
11

Figure 9: Spectral Plot of Structural Brain Net-Figure 10: Spectral Plot of Structural Brain Network.(Age: 9, Region 80, λ2 = 0.2013)
work.(Age: 9, Region 500, λ2 = 0.1323)

Figure 11: Spectral Plot of Structural Brain Network.Figure 12: Spectral Plot of Structural Brain Network.
(Age: 67, Region 80, λ2 = 0.2662)
(Age: 67, Region 500, λ2 = 0.1506)

Using the parameters introduced above, we explore
their statistical variations for the two comparison groups
as a function of the number of nodes in the network. Table
6 shows the p-values of the samples of the spectral parameters corresponding to the two sets of population groups
defined in Table 3. It is easy to see that the modularity
parameter and second eigenvalue are significantly different
for each set of population groups, especially for large-scale
structural brain networks. These findings suggest that
fine-grained brain networks are better able at capturing
subtle structural patterns than coarser-grain networks.

spectrum seems to be more robust to the threshold value
chosen to control the sparsity of brain networks. Table
8 lists the values of some graph-theoretic measures of a
500-node brain network using different threshold values.
As can be seen from these results, the values of the graph
theoretic parameters change substantially as a function of
the threshold value, which cast some doubt on their utility
to analyze brain networks.
Table 8: Robustness of Graph Measurements with Different Threshold Value


0
0.0001
0.001
0.01

5.2.3. Robustness of Spectral Graph Parameters
In this section we shed some light on the relative robustness of the two types of network analysis. It is wellknown that the graph theoretic measures depend heavily on the preprocessing step. Recall that the preprocessing step eliminates edges with very weak connections
and transforms relatively dense graphs into sparse graphs.
Graph measures such as the characteristic path length and
clustering coefficient are sensitive to the sparsity ratio defined by the fraction of non-zero entries in the adjacent
matrix, which is controlled by the threshold value used
during the preprocessing step. However, the Laplacian

CPL
1.2819
1.4128
2.0083
3.4294

Eglobal
0.8590
0.7936
0.5499
0.3328

C
0.8185
0.7075
0.4294
0.4209

Sparsity ratio
0.7181
0.5873
0.1538
0.0333

Figure 13 shows the spectral plots corresponding to the
weighted brain networks using the same set of threshold
values as before. The spectral plots are similar despite that
the fact that the sparsity ratios are significantly different.
These results indicate that the spectral properties of the
high-resolution brain networks are more robust than tra12

Figure 14: Time costs of brain parcellation scheme for one iteration

of an arbitrary subject, we obtain a connectivity graph
characterized in Table 9.
Figure 5.3 shows the time costs (in seconds) of each major step for a single iteration of our parcellation algorithm.
Note that our algorithm typically converges in less than 4
iterations. Notice that the construction of the similarity matrix and the execution of the k-means algorithm are
substantially accelerated, which leads to about 5x speedup
of the overall execution time compared to the fast Matlab implementation. We also obtain a good speedup for
computing the k eigenvalues and eigenvectors making use
of very fast sparse matrix-vector multiplication. On our
CPU-GPU system, the overall execution time of each iteration is reduced to less than 10 minutes for k = 500, which
makes it quite feasible to generate almost arbitrarily large
size parcellations.

Figure 13: Spectral Plot of Brain Network (Regions: 500), Subject
Age: 41, Gender: male

Table 9: Problem size

# of voxels
142541

# of edges
3992290

# of clusters
500

ditional graph theoretic measures relative to the threshold
values using during preprocessing.
5.2.4. Classification based on Network Parameters
Tables 5 and 6 illustrated the increasingly statistical
differences of graph measures for our distinct population
groups as we increase the number of nodes. We provide
here an alternative evidence by building classifiers for each
of our pair of population groups and measure their accuracy relative to the number of nodes in the network. We
train a Support Vector Machine (SVM) classifier with linear kernel7 on the two sets of subjects of Comparison 1
and Comparison 2 defined in table 3 using a variety of
graph features. To measure the classification accuracy, we
use the Leave-one-out cross-validation (LOOCV) on both
sets of population groups. Table 7 shows the classification
accuracy based onthe use of different graph theoretic and
graph spectral features.
For the fine-grained brain networks (k = 500), the SVM
classifier achieve consistently better accuracy than using
coarse-grain structural brain networks. Note that the relative high accuracy achieved by using only the second eigenvalue of the normalized Laplacian of the network.

6. Conclusion
In this paper, we presented a parallel connectivitybased whole-brain parcellation scheme that was implemented
on a heterogeneous CPU-GPU platform. The resulting
algorithm generates large-scale parcellations quickly with
provably structurally homogeneous regions. The algorithm
was used to build structural brain networks with varying
numbers of nodes up to 500 nodes. These brain networks
were then analyzed using traditional graph theoretic measures as well as new measures that are based on spectral graph theory. Using data from populatioin groups
that are expected to have significantly different structural
patterns, we illustrated the increasing statistical power of
the neworks as the number of nodes grows. Moreover, we
provided preliminary evidence on the robustness of using
spectral parameters to analyze brain networks and their
potential for discriminating between different population
groups. We are currently persuing this line of research
trying to shed more light on the use of spectral graph theory to characterize brain networks for different population
groups.

5.3. Computational Efficiency
In this section, we briefly provide an indication of the
time efficiency of our parallel implementation of the connectivity based brain parcellation scheme. We compare
our fast implementation to the best multithreaded implementations we could achieve using Matlab and Python scientific packages. After removing isolated voxels of the data
7 http://www.mathworks.com/help/stats/svmtrain.html

13

7. Acknowledgment
[17]
We gratefully acknowledge funding provided by The
University of Maryland/Mpowering the State through the
Center for Health-related Informatics and Bioimaging (CHIB)
and the NVIDIA Research Excellence Center at the Uni[18]
versity of Maryland.

References
[19]

[1] E. Bullmore, O. Sporns, Complex brain networks: graph theoretical analysis of structural and functional systems, Nature
Reviews Neuroscience 10 (3) (2009) 186–198.
[2] O. Sporns, Structure and function of complex brain networks,
Dialogues in clinical neuroscience 15 (3) (2013) 247.
[3] Q. Wang, R. Chen, J. JaJa, Y Jin, L. Hong, and E. Herzkovits,
“Connectivity-Based Brain Parcellation: A Connectivity-Based
Atlas for Schizophrenia Research”, submitted to Neuroinformatics.
[4] E. H. Herskovits, L. E. Hong, P. Kochunov, H. Sampath,
R. Chen, Edge-centered dti connectivity analysis: Application
to schizophrenia, Neuroinformatics (2015) 1–9.
[5] M. Kaiser, A tutorial in connectome analysis: topological and
spatial features of brain networks, Neuroimage 57 (3) (2011)
892–907.
[6] A. Zalesky, A. Fornito, I. H. Harding, L. Cocchi, M. Yücel,
C. Pantelis, E. T. Bullmore, Whole-brain anatomical networks:
does the choice of nodes matter?, Neuroimage 50 (3) (2010)
970–983.
[7] M. Ingalhalikar, A. Smith, D. Parker, T. D. Satterthwaite, M. A.
Elliott, K. Ruparel, H. Hakonarson, R. E. Gur, R. C. Gur,
R. Verma, Sex differences in the structural connectome of the
human brain, Proceedings of the National Academy of Sciences
111 (2) (2014) 823–828.
[8] G. Gong, P. Rosa-Neto, F. Carbonell, Z. J. Chen, Y. He,
A. C. Evans, Age-and gender-related differences in the cortical anatomical network, The Journal of neuroscience 29 (50)
(2009) 15684–15693.
[9] Y. Iturria-Medina, Anatomical brain networks on the prediction
of abnormal brain states, Brain connectivity 3 (1) (2013) 1–21.
[10] M. Rubinov, O. Sporns, Complex network measures of brain
connectivity: uses and interpretations, Neuroimage 52 (3)
(2010) 1059–1069.
[11] N. Tzourio-Mazoyer, B. Landeau, D. Papathanassiou, F. Crivello, O. Etard, N. Delcroix, B. Mazoyer, M. Joliot, Automated
anatomical labeling of activations in spm using a macroscopic
anatomical parcellation of the mni mri single-subject brain,
Neuroimage 15 (1) (2002) 273–289.
[12] Q. Cao, N. Shu, L. An, P. Wang, L. Sun, M.-R. Xia, J.-H.
Wang, G.-L. Gong, Y.-F. Zang, Y.-F. Wang, et al., Probabilistic
diffusion tractography and graph theory analysis reveal abnormal white matter structural connectivity networks in drug-naive
boys with attention deficit/hyperactivity disorder, The Journal
of Neuroscience 33 (26) (2013) 10676–10687.
[13] S. B. Eickhoff, K. E. Stephan, H. Mohlberg, C. Grefkes, G. R.
Fink, K. Amunts, K. Zilles, A new spm toolbox for combining probabilistic cytoarchitectonic maps and functional imaging
data, Neuroimage 25 (4) (2005) 1325–1335.
[14] D. Moreno-Dominguez, A. Anwander, T. R. Knösche, A hierarchical method for whole-brain connectivity-based parcellation,
Human brain mapping 35 (10) (2014) 5000–5025.
[15] M. Ruschel, T. R. Knösche, A. D. Friederici, R. Turner,
S. Geyer, A. Anwander, Connectivity architecture and subdivision of the human inferior parietal cortex revealed by diffusion
mri, Cerebral Cortex 24 (9) (2014) 2436–2448.
[16] R. B. Mars, S. Jbabdi, J. Sallet, J. X. O’Reilly, P. L. Croxson,
E. Olivier, M. P. Noonan, C. Bergmann, A. S. Mitchell, M. G.
Baxter, et al., Diffusion-weighted imaging tractography-based
parcellation of the human parietal cortex and comparison with

[20]

[21]

[22]

[23]

[24]

[25]

[26]
[27]

[28]

[29]

[30]

[31]

[32]

[33]

[34]

[35]

[36]

14

human and macaque resting-state functional connectivity, The
Journal of Neuroscience 31 (11) (2011) 4087–4100.
R. B. Mars, J. Sallet, U. Schüffelgen, S. Jbabdi, I. Toni, M. F.
Rushworth, Connectivity-based subdivisions of the human right
temporoparietal junction area: evidence for different areas participating in different cortical networks, Cerebral cortex 22 (8)
(2012) 1894–1903.
J. Sallet, R. B. Mars, M. P. Noonan, F.-X. Neubert, S. Jbabdi,
J. X. O’Reilly, N. Filippini, A. G. Thomas, M. F. Rushworth, The organization of dorsal frontal cortex in humans and
macaques, The Journal of Neuroscience 33 (30) (2013) 12255–
12274.
F.-X. Neubert, R. B. Mars, A. G. Thomas, J. Sallet, M. F.
Rushworth, Comparison of human ventral frontal cortex areas
for cognitive control and language with areas in monkey frontal
cortex, Neuron 81 (3) (2014) 700–713.
F.-X. Neubert, R. B. Mars, J. Sallet, M. F. Rushworth, Connectivity reveals relationship of brain areas for reward-guided
learning and decision making in human and monkey frontal cortex, Proceedings of the National Academy of Sciences 112 (20)
(2015) E2695–E2704.
A. Anwander, M. Tittgemeyer, D. Y. von Cramon, A. D.
Friederici, T. R. Knösche, Connectivity-based parcellation of
broca’s area, Cerebral Cortex 17 (4) (2007) 816–825.
S. Mittal, J. S. Vetter, A survey of cpu-gpu heterogeneous computing techniques, ACM Computing Surveys (CSUR) 47 (4)
(2015) 69.
C. Lee, W. W. Ro, J.-L. Gaudiot, Boosting cuda applications
with cpu–gpu hybrid computing, International Journal of Parallel Programming 42 (2) (2014) 384–404.
J. Agulleiro, F. Vazquez, E. Garzon, J. Fernandez, Hybrid computing: Cpu+ gpu co-processing and its application to tomographic reconstruction, Ultramicroscopy 115 (2012) 109–114.
J. Wu, J. JaJa, Optimized fft computations on heterogeneous
platforms with application to the poisson equation, Journal of
Parallel and Distributed Computing 74 (8) (2014) 2745–2756.
J. Wu, J. JaJa, Achieving native gpu performance for out-ofcard large matrix multiplication, Parallel Processing Letters.
K. Li, W. Yang, K. Li, A hybrid parallel solving algorithm
on gpu for quasi-tridiagonal system of linear equations, IEEE
Transactions on Parallel and Distributed Systems.
E. N. Akimova, D. V. Belousov, Parallel algorithms for solving
linear systems with block-tridiagonal matrices on multi-core cpu
with gpu, Journal of Computational Science 3 (6) (2012) 445–
449.
Y. Jin, J. JaJa, A high performance implementation of spectral clustering on cpu-gpu platforms, in: International Parallel
and Distributed Processing Symposium Workshop (IPDPSW),
2016.
R. B. Lehoucq, D. C. Sorensen, C. Yang, ARPACK users’
guide: solution of large-scale eigenvalue problems with implicitly restarted Arnoldi methods, Vol. 6, SIAM, 1998.
R. B. Lehoucq, D. C. Sorensen, Deflation techniques for an
implicitly restarted arnoldi iteration, SIAM Journal on Matrix
Analysis and Applications 17 (4) (1996) 789–821.
D. C. Sorensen, Implicit application of polynomial filters in akstep arnoldi method, Siam journal on matrix analysis and applications 13 (1) (1992) 357–385.
Y. Jin, J. F. JaJa, R. Chen, E. H. Herskovits, A data-driven approach to extract connectivity structures from diffusion tensor
imaging data, in: 2015 IEEE International Conference on Big
Data, IEEE, 2015.
J. L. Andersson, M. Jenkinson, S. Smith, et al., Non-linear registration, aka spatial normalisation fmrib technical report tr07ja2,
FMRIB Analysis Group of the University of Oxford.
T. Behrens, H. J. Berg, S. Jbabdi, M. Rushworth, M. Woolrich,
Probabilistic diffusion tractography with multiple fibre orientations: What can we gain?, Neuroimage 34 (1) (2007) 144–155.
T. Behrens, M. Woolrich, M. Jenkinson, H. Johansen-Berg,
R. Nunes, S. Clare, P. Matthews, J. Brady, S. Smith, Characterization and propagation of uncertainty in diffusion-weighted mr

[37]
[38]

[39]

[40]

[41]

[42]
[43]
[44]

[45]

imaging, Magnetic resonance in medicine 50 (5) (2003) 1077–
1088.
U. Von Luxburg, A tutorial on spectral clustering, Statistics
and computing 17 (4) (2007) 395–416.
A. Y. Ng, M. I. Jordan, Y. Weiss, et al., On spectral clustering: Analysis and an algorithm, Advances in neural information
processing systems 2 (2002) 849–856.
S. X. Yu, J. Shi, Multiclass spectral clustering, in: Computer
Vision, 2003. Proceedings. Ninth IEEE International Conference on, IEEE, 2003, pp. 313–319.
Y. Sun, R. Lee, Y. Chen, S. Collinson, N. Thakor, A. Bezerianos, K. Sim, Progressive gender differences of structural
brain networks in healthy adults: a longitudinal, diffusion tensor imaging study, PloS one 10 (3) (2015) e0118857.
Q. Wang, R. Chen, J. JaJa, Y. Jin, L. E. Hong, E. H. Herskovits,
Connectivity-based brain parcellation, Neuroinformatics 14 (1)
(2016) 83–97.
V. Latora, M. Marchiori, Efficient behavior of small-world networks, Physical review letters 87 (19) (2001) 198701.
F. R. Chung, Spectral graph theory, Vol. 92, American Mathematical Soc., 1997.
A. Banerjee, The spectrum of the graph laplacian as a tool
for analyzing structure and evolution of networks, Ph.D. thesis
(2008).
S. C. de Lange, M. A. de Reus, M. P. van den Heuvel, The
laplacian spectrum of neural networks, Frontiers in computational neuroscience 7 (1).

15

