To allow for more straightforward estimations of outlierness, as well as similarities and differences between cortical vertices, we applied an invertible non-linear transformation to map the features into a latent representation where the PDF follows an approximately multivariate normal distribution using the Rotation-Based Iterative Gaussianization (RBIG) procedure \cite{Laparra_2011}.  RBIG is a fast iterative procedure consisting of a sequence of pairs of transformations: 1) a non-linear transformation applied to each of the columns (marginals) of the data matrix and 2) a linear transformation applied to the entire data matrix.  The non-linear column-wise operation is a univariate Gaussianization that converts percentile scores computed using the rank transformation to standard scores (scikit-learn's  QuantileTransformer). The orthogonal transformation is performed using PCA, with all components retained after each iteration. Ten iterations were performed. 

Feature Comparison

Because this model theoretically represents an efficient and overcomplete representation of the local features of each image, we hypothesized that we should be able to predict the values of other, more commonly used local features such as curvature, sulcal depth, cortical thickness, and gray/white contrast (as calculated using FreeSurfer) or measures of myelination, here calculated by dividing the Tintensity by the T2 intensity, as in  \cite{Glasser_2011}, then sampled onto the gray-white function surface.  Using our feature set in as the input, quadratic regression models (ordinary least-squares) were estimated using scikit-learn for each feature.  using a 10-fold cross validation procedure.  Each fold was trained on all cortical vertices from 25 randomly selected healthy volunteers and tested on all cortical vertices from the remaining 5 subjects.  Performance was evaluated for each model using the coefficient of determination \(r{^2}\); effect size was reported as in \cite{cohen1988statistical}\(r=0.1\) as small\(r=0.3\) as medium, and \(r=0.5\) as large.

Distance and Similarity Estimation

To smooth the data, we created patches centered at every vertex, including all neighboring vertices within a 5 mm radius.  Across subjects, homotopic and heterotopic patches were defined as patches in the same or different locations, respectively. For each patch \(p\), the center \(\mathbf{m}_p\) was computed by averaging the feature vectors of the vertices within the patch: \(\mathbf{m}_p = \sum_{\mathbf{x}\in p} \mathbf{f}(\mathbf{x})\).  The direction \(\mathbf{\hat{u}}_p\), the unit vector pointing in the direction of the patch's center, was computed by dividing the patch's mean feature vector by its length \(\hat{u}_p = \mathbf{m}_p / \| \mathbf{m}_p \|\).  The probability of finding a patch with a given average feature vector \(\mathbf{m}_p\) depends only on the magnitude of the feature vector \(\| \mathbf{m}_p \|\), which is also the Mahalanobis distance, \(d\), defined as the distance from the origin to the patch's center with \(d^2\) following a cumulative chi-squared distribution.  The similarity between two patches \(p_i\) and \(p_j\) can be assessed using simple metrics, such as 1) the Euclidean distance between their centers \(d_{ij} = \| \mathbf{m}_i - \mathbf{m}_j \|\), and 2) the cosine similarity between their directions \(s_{ij} = \mathbf{\hat{u}_i} \cdot \mathbf{\hat{u}_j}\).

Global Anomaly Detection

We hypothesized that cortical lesions, but also possibly some normally "atypical" cortical regions, would appear as global outliers in our feature space. To identify such regions of normal cortex, we calculated the average Mahalanobis distance for each cortical patch across all HVs (fig).  Specific outlier regions were identified by thresholding the average distance map across HVs at a threshold of 2.7, retaining 4.3% of the patches ( \(p=0.043\)), and clustering with a minimum of 30 nodes.  As exemplars, for further analysis we selected 2 of the resulting outlier regions of interest (ROI) in the anterior insula and primary motor cortex.

Similarity Estimation

Although some brain regions appear to be consistent global outliers across HVs, this does not mean that these regions are similar to each other.  In our representation, similarities or differences in combinations of features can be represented as differences in direction.  To explore this in our 2 exemplar outlier ROIs, we defined the center and direction of each ROI as the average of the patches within that ROI across all HVs.  We similarly defined the average FCD ROI center and direction by averaging all of the patches within each MRI+ FCD mask (n=10), then averaging across the FCDs, to give each lesion equal weighting. Using the average direction of the insula ROI as the x-axis and the average direction of the motor ROI as the y-axis, we plotted the patches within each outlier ROI, the FCD ROI center, and 1000 randomly selected cortical patches for comparison (fig).  
To provide an overall measure of the 14-dimensional similarity of all cortical patches to each selected ROI, we computed the scalar projection of the mean feature vector for each patch onto the ROI average unit vector.  For our 3 ROIs, we calculated these similarity maps for each healthy control and displayed the average similarity map for each ROI on the cortical surface (fig).

FCDs as Outliers

We next wished to assess to what extent FCD lesions in MRI+ patients were global outliers in our feature space. We compared the average Mahalanobis distances of 1000 randomly selected patches from 1) normal HV cortex, 2) insular and motor cortex outlier ROIs, and 3) FCD lesion masks.  As FCDs were not more anomalous and appear similar to some portions of normal cortex, , we then hypothesized that to be identifiable, FCDs must differ from the expected appearance of their homotopic region. We implemented a local normalization procedure by removing the local mean of homotopic patches calculated across all HVs, followed by variance normalization, with the local standard deviation for each patch calculated across HVs after removal of the local mean for each control subject individually.  We then again compared the Mahalanobis distances of normal cortex, the outlier ROIs, and FCDs.
Although we expected that local normalization should minimize reproducible regional differences across normal cortex, the effect on FCD similarity (or dissimilarity) to normal cortex and within the FCDs as a class was unknown.  Therefore, before and after local normalization, we computed the Euclidean distance and cosine similarity between 1000 randomly selected FCD patches and 1) 1000 homotopic patches (within the FCD mask) in HVs, 2) 1000 heterotopic patches in HVs, and 3) 1000 other non-overlapping FCD patches across patients.  We compared these distances to those between 1000 pairs of homotopic patches and 1000 pairs of randomly selected heterotopic patches across HVs.  For each comparison, we report the mean and standard deviation;  Welch's t-test was used to compare the distribution of distances and cosine similarities between patches, as they were approximately normally distributed but often had unequal variances.  Effect size was reported as in \cite{cohen1988statistical}\(d=0.2\) as small\(d=0.5\) as medium, and \(d=0.5\) as large.
Finally, to visualize the variability across FCD lesions before and after local normalization, we used a singular value decomposition (SVD) to select the 2 components explaining the most variance across the FCDs and used these as the x- and y-axes to create a 2D plot of the location of each FCD patch, along with the center of each individual FCD and the average FCD center across all subjects.  We also plotted the centers of 1000 randomly sampled healthy control patches for comparison.

FCD Detection

Finally, we created an automated FCD detection algorithm. As in \cite{Hong2014}, we used a two stage method to improve specificity, using a quadratic discriminant analysis implemented in Scikit-learn, using a leave-one-out cross-validation strategy for each FCD patient and each healthy volunteer to build a subject specific model trained using the labeled data from the other subjects as follows:
  1. For the first stage, we randomly sampled with replacement: a) 1000 normal patches evenly across the healthy volunteers, and b) 1000 FCD patches from the patient lesions masks, and fit a quadratic discriminant analysis (QDA) model to classify patches as FCD or HV.
  2. We applied this model to all of the patches in all of the training subjects and retained all clusters of greater than 5 patches classified as FCD.
  3. For the second stage, we sampled with replacement from clusters surviving the first stage:  a) 1000 HV patches and b) 1000 FCD patches, and fit a second QDA model.
  4. We applied this second model to all of the patches that survived the first stage model in each test subject, retaining all clusters of greater than 5 patches classified as FCD.
  5. For each surviving cluster, we computed a cluster weight, calculated as area \(\times\)mean FCD probability, as well as the cluster's rank among surviving clusters in that subject.
Sensitivity of the classifier was assessed initially for MRI+ patients and then for all patients.  We calculated a receiver operating characteristics (ROC) curve and area under the curve (AUC).  True positives (or detections) were defined as MRI+ patients with co-localization of a detected cluster and the manual lesion mask, and for MRI- patients, as overlap of a detected cluster with the resection mask.  False positives were defined as HVs with detected clusters detected above a given threshold.  The optimal threshold to apply to the resulting clusters was determined using the Youden Index (calculated as sensitivity + specificity - 1).  Additionally, lesion detection and number of extralesional clusters (outside of the lesion or resection masks or in HVs) were assessed for all patients and HVs individually.  

Results

Study Participants

A total of 15 patients with drug resistant focal epilepsy and FCD (median age 27, range 15-53, 11 females) and 30 healthy controls (median age 23, range 8-63, 12 females) were included in this study.  FCDs were identified in the radiological reports in 6/15 patients, an additional 5 patients had lesions identifiable on post-hoc analysis (MRI+ n=11).  Patient demographics and further details are reported in Table (summary table).

Feature Comparison

Using our novel feature set, we were able to predict several several features commonly used to describe the cortical variability .  Our prediction of the Freesurfer measures of curvature and sulcal depth were highly accurate (r= 0.81 and r= 0.72 respectively, both large effect sizes).  This was expected, given that our Laplacian filters are a measure of local curvature, and sulcal depth is highly correlated with curvature.  Also accurately predicted with large effect sizes were gray-white contrast (r2 = 0.65), represented in our feature set in the form of the gradient filters, myelin contrast (r2 = 0.54), derived from T1 and T2 intensities, and cortical thickness (r2 = 0.56), which is less directly represented in our initial feature set.

Global Anomaly Detection and Similarity Estimation

We identified 11 and 8 consistently outlying clusters in the left and right hemispheres respectively across HVs, located in portions of the primary sensorimotor cortices, as well as limbic and paralimbic regions (fig).  Patches from two selected outlier regions in insular and primary motor cortex can be seen to lie in relatively well separated clusters centered around the average ROI center, distinct both from each other and from normal cortex, with the average FCD appearing quite close to the average insula vector.  Similarity maps demonstrate that these differences and similarities carry into our 14 dimensional feature space, with insular and FCD ROI similarity maps appearing quite similar,  clearly differing from the spatial distribution of the motor cortex ROI similarity map.

FCDs as Outliers

 We next wished to evaluate whether FCDs are global outliers in our feature space compared to normal cortex, particularly those areas previously identified as outlier regions.  FCD patches do in fact appear to be outliers compared to randomly sampled "normal" cortical patches with a large effect size ($3.24\pm0.64$ vs $2.61\pm0.57$, $d=1.041$), but not more so than other normal outlying ROIs (motor $3.45\pm0.43$, $d=0.40$, insula $3.37\pm0.36$, $d=0.25$).  Following local normalization, as expected, patches in the "normal" outlying ROIs on average moved closer to the center of the distribution, appearing more similar to randomly selected normal cortical patches than previously (normal: $0.58\pm0.15$ versus motor: $0.66\pm0.17$ ($d=0.48$), and insula: $0.59\pm0.17$, $d=0.10$), as these regions are fairly reproducibly spatially located in healthy volunteers and thus now appear more "normal" for their location.  In contrast, following local normalization, FCDs remain as significant outliers ($0.93\pm0.25$) compared to both normal cortex and the outlying "normal" ROIs with large effect sizes (vs normal cortex $d=1.72$, motor $d=1.32$, insula $d=1.62$). 
Using our normative model, we explored several aspects of normal cortical variability.  First, we calculated our similarity metrics of Euclidean distance ($d_e$) and cosine similarity ($cos$) in HVs between: 1) pairs of homotopic patches across subjects  ($d_e$=$2.67\pm0.73$, $cos$=$0.45\pm0.28$), 2) heterotopic patches across subjects ($d_e$=$3.69\pm0.79$, $cos$=$0.00\pm0.29$).  In both distance and cosine similarity, homotopic patches appear more similar to each other than heterotopic patches (p<0.001, $d=1.35$ and $d=1.59$ for $d_e$ and $cos$, respectively) with a large effect size.  As expected, these differences between homotopic and hetereotopic patches are less notable after carrying out local normalization, with $d_e$ = $2.67\pm0.73$ vs. $d_e$=$2.65\pm0.65$ ($d=0.04$), and $cos$=$-0.03\pm0.34$ vs. $cos$=$0.00\pm0.29$ ($d=0.11$). 
Assessing the similarity of FCD patches to normal cortex, prior to local normalization....
following local normalization the average distance between pairs of FCD patches and normal cortex appears similar regardless of the location of the cortex ($d_e$=$3.63\pm0.94$ and$d_e$=$3.61\pm0.92$, between pairs of FCD patches and homotopic or heterotopic patches, respectively ($d=0.02$)). Interestingly, these distances are quite similar to the within class distance of FCD patches across individuals ($3.48\pm0.96$, $d=0.16$ and $d=0.14$ for FCD vs homotopic and heterotopic patches, respectively).  For cosine similarity, there is again little difference when comparing FCD to homotopic ($cos$=$0.00\pm0.30$) or heterotopic ($cos$=$0.00\pm0.28$) patches, but cosine similarity across FCDs ($0.39\pm0.26$) is much greater than between FCDs and normal cortex ($d=1.38$ and $d=1.41$ for FCD vs homotopic and heterotopic patches, respectively).  This suggests that FCDs are located in a relatively consistent direction, but are not clustered tightly, as their distances from each other are similar to that of randomly selected pairs of patches.  This can be observed in graphic form in figure XX, in which the center of each FCD is shown, as well as the average center and direction of the average FCD vector.  Additionally, in figure XX, the MRI appearance within individual FCDs are demonstrated,  as are the locations of the lesional patch centers in 2D, highlighting the variable appearances within and across FCDs both visually and in this 2D space.