Geometric Deep Learning on Multiple Spectral Embeddings: Learning of Cortical Surface Data

Let brain graph $\mathcal{G}=\{\mathcal{V},\mathcal{E}\}$ from set of $N$ vertices with position $\mathbf{x}_{i}\in\mathbb{R}^{3}$, and edges of a surface model $S$. The normalized graph Laplacian operator with weighted adjacency matrix $\mathbf{A}$ and diagonal degree matrix $\mathbf{D}$ is of the form, $\mathbf{L}=\mathbf{I}-\mathbf{D}^{-1}\mathbf{A}$. Then, the spectral decomposition of this graph is $\mathbf{L}=\mathbf{U}\mathbf{\Lambda}\mathbf{U}^{-1}$. The vibration modes of the surface model $S$ is depicted by the sorted eigenvalues $\mathbf{\Lambda}=\mathrm{\textbf{Diag}}(\lambda_{0},\ldots,\lambda_{N})$ and their associated eigenfunctions $\mathbf{U}=[\mathbf{u}_{0},\ldots,\mathbf{u}_{N}]$. The n-dimensional spectral coordinates in normalized scale is $\widehat{\mathbf{U}}=\mathbf{\Lambda}^{\frac{1}{2}}\mathbf{U}$. Unlike Euclidean domain for navigating on surfaces, spectral coordinates can aid in moving over the surface for learning.

The Laplacian for different brain graph differs due to the inconsistency in adjacency matrix across subjects. Hence, the spectral basis for multiple brains graphs will be misaligned and the learning algorithm would be incapable of generalization. We handle this issue by realigning all spectral representations to an arbitrary reference. We can compute a spectral transformation $T^{(f\rightarrow g)}$ that can produce a new spectral basis $\mathbf{U}_{T}(f)$ for a brain graph $\mathbf{U}(f)$ matching to a reference basis $\mathbf{U}(g)$ such that $T^{(f\rightarrow g)}=(\mathbf{U}(g))^{-1}\mathbf{U}_{T}(f)$. The first five spectral coordinates with Iterative closest point, in practice, are sufficient for matching spectral representation of two different brains. The aligned spectral coordinates $spe(n)$ thus obtained for a vertex x for a brain $g$ is $n^{th}$ row of the matrix $\widehat{\mathbf{U}}T^{(g\rightarrow ref)}$.

We propose a deep learning framework to learn on graphs by adapting the traditional convolution principles applied on rigid Euclidean space. A point $x$ on a cartesian grid has a defined neighborhood at a fixed distance. Consider nodes $i$ and $j$ from the built graph $\mathcal{G}=\{\mathcal{V},\mathcal{E}\}$. from set of vertices with position x, and edges of a surface model $S$. The geodesic distance between vertices of a surface model vary depending on the folds and valleys. A fixed grid based neighborhood cannot be defined for a graph structure due to its inherent properties. Hence, we propose to use a learnable gaussian kernel to measure the influence of neighboring nodes on a given node $i$. A Gaussian graph kernel between node $i$ and node $j$ is defined as,

where $g_{ij}$ indicates node connectivity of $i$ and $j$; $X_{i}$ is the spectral embedding of the node $i$; $\mu_{k}$ and $\sigma_{k}$ are the Gaussian parameters of the filter kernel learned over epochs.

A graph convolution filter $w_{pqk}^{(l)}$ between input $q$ and output $p$ is learned producing a $M_{l}$ filter response maps at a layer $l$. The graph convolution operation is defined as,

where, $b_{q}^{(l)}$ is the bias for kernel $k$ for input feature $q$; $y_{jq}^{(l)}$ is the input feature $q$ of node $j$.

The proposed architecture contains three graph convolution layers ($l=3$) with feature maps $M_{1}=32;M_{2}=64;M_{3}=32$. The size of the last layer corresponds to the number of parcels to be segmented (32 in this case). Each layer has four Gaussian kernels $K=4$ operating on $N$ connected neighbours of node $i$. A non linear activation function, for instance, Leaky ReLU is applied after each layer to obtain filter response $y_{ip}^{(l)}=\max(0.01z_{ip}^{(l)},\,z_{ip}^{(l)})$.

Since the parcels to be segmented are mutually exclusive, we choose softmax layer at end of last graph convolution layer to obtain a probability of belonging to a parcel for all node. The softmax function is given by $\frac{exp(z_{ip}^{l})}{\sum_{q}exp(z_{iq}^{l})}$. The error $E$ between the predicted $y_{ic}$ and the target $q_{ic}$ is measured using the cross-entropy loss given by

This loss is minimized by propagating the error using standard gradient descent optimization scheme. Along with graph convolution filter weights, the parameter $\mu_{k}$ and $\sigma_{k}$ of the Gaussian kernel is updated.

This experiment evaluates the improvement of moving the learning operations from the Euclidean domain to a Spectral domain. To do so, we compare the classification accuracy on 32 cortical parcels when running our algorithm, respectively, in the Euclidean and Spectral domains. Quantitative results are measured in terms of average Dice overlap ($2|A\cap B|/(|A|+|B|)$) and Hausdorff distances. Qualitative results are shown in Fig. 2.

Euclidean domain – In our baseline, similarly to the latest approaches of graph convolutions networks [Monti2017Geometric], we learn from input features in the Euclidean domain. Each cortical point is represented using sulcal depth and its spatial location, i.e., $f_{i}^{eu}=(S_{d}^{i},X,Y,Z)$. The architecture remains the same as described earlier (Fig. 1). We train our network on the randomly split training set until the performance decreases on the validation set. The average Dice overlap across all parcels in our dataset is 59% ($\pm$ 20.5, min/max = 0/85.54%). The Hausdorff distance averaged across all parcels is 18.87mm. Fig. 2 clearly illustrates the current limitation of existing graph convolutions approaches. Spatial ambiguities are inherent in the Euclidean domain, i.e., neighboring coordinates in space may not necessarily be neighbors on surfaces. This ambiguity may confuse training on highly convoluted surfaces such as the brain, and possibly explains the strong spatial irregularities observed in the parcels boundaries.

Spectral domain – Our contribution is to operate in a geometry-aware Spectral domain. We now learn on input features represented in the Spectral domain, where each cortical point consists of the sulcal depth with the first three spectral coordinates, i.e, $f_{i}^{sp}=(S_{d}^{i},spe(1),spe(2),spe(3))$. We use the same architecture and data split as before. The average Dice overlap across all parcels improves to 82% ($\pm$ 5.32, min/max = 69.45/95.09%). Details for each parcels are shown in Fig. 3. The Hausdorff distance averaged across all parcels is now reduced to 12.30mm. This is a 40% improvement over learning in the conventional Euclidean domain. The qualitative results of Fig. 2 show that our cortical parcellation is almost similar to the manual parcellation. The boundary, however, is irregular and requires further regularization.

As an illustration of further refinement, we apply a standard graphcut [graphcut] regularization for both, Euclidean and Spectral outputs. To do so, the graphcut algorithm is set using the surface adjacency matrix as graph edge weights, and the output node probabilities as graph node weights. The graphcut regularization further improves the overall classification accuracy from 68.24% to 71.08% in the Euclidean domain, and from 84.35% to 85.35% in the Spectral domain. Similar improvement is observed in terms of Hausdorff distance, with a reduction from 18.87mm to 15.96 in the Euclidean domain, and from 12.30mm to 10.62mm in the Spectral domain. The Dice overlap scores remain, however, similar after regularization, perhaps due to a correct initial overlap of parcels by our algorithm.

In this experiment, we illustrate the benefits of exploiting neighborhoods of surface data over pointwise information in a spectral domain. To follow standard approaches of learning pointwise surface data [Lombaert2015Spectral], we set our weighted adjacency matrix as an identity matrix. This removes any neighboring connection. The graph convolutions, therefore, requires no update for $\sigma$ and $\mu$. Training an pointwise information, $f_{i}^{sp}=(S_{d}^{i},spe(1),spe(2),spe(3))$, results in an average Dice overlap of 75% ($\pm$ 14.3%, min/max = 16.16/92.32%). As a reminder, training on neighborhoods results in 82% of Dice score. Exploiting neighborhoods of surface data improves, therefore, the performance of our algorithm. A closer look to the performance scores for each parcels (Fig. 3) also reveals a general improvement when exploiting neighborhoods over pointwise surface data. Results for all parcels are also at par with FreeSurfer’s performance.

In the MindBoggle dataset, the manual parcellations were created by experts who corrected parcels initially generated from FreeSurfer. The agreement between manually corrected parcels and FreeSurfer is $83\%$ ^{fancyline,backgroundcolor=red!25,bordercolor=redfancyline,backgroundcolor=red!25,bordercolor=red}todo: fancyline,backgroundcolor=red!25,bordercolor=redHervé says: standard deviation? of Dice score. Our algorithm produces a 82% Dice score with respect to manual parcellation, while gaining a significant improvement in computation time, from hours to seconds for processing one subject. ^{fancyline,backgroundcolor=red!25,bordercolor=redfancyline,backgroundcolor=red!25,bordercolor=red}todo: fancyline,backgroundcolor=red!25,bordercolor=redHervé says: We’re missing a note on computation time here Our graph convolutional framework also results in a stable performance when compared to an Euclidean approach or to FreeSurfer, who respectively have a standard deviation on Dice scores of 20.5% and 10.3% across all parcels. Our spectral approach varies by 5.3%.