Multimodal Tissue Segmentation is better

doi:10.21203/rs.3.rs-4870185/v1

Download PDF

Research Article

Multimodal Tissue Segmentation is better

https://doi.org/10.21203/rs.3.rs-4870185/v1

This work is licensed under a CC BY 4.0 License

You are reading this latest preprint version

Being able to distinguish tissues and quantify features in an image is ubiquitous in medical imaging, allowing, for instance, the study of differences between clinical groups or the investigation of the impact of an intervention. Assigning the tissue type is also a fundamental preprocessing step in many neuroimaging applications like image registration, normalisation or even simple masking. Because of such ubiquity, thorough investigations of segmentation algorithms are necessary to determine in which conditions they work best. The SPM implementation of tissue segmentation is a commonly used tool in this context, providing voxel-wise probabilistic estimates of brain grey and white matter tissues, cerebrospinal fluid, soft tissues, and bones. Different estimates of tissue density and/or volumes have, however, been observed using unimodal vs. multimodal inputs. Here, we contend that possible misinterpretations arise from mis-specifying parameters of the generative model underlying tissue segmentation. Using T1 weighted vs. T1 and T2 weighted images as input while also varying the number of Gaussians (1 vs. 2 for brain tissues) used in the generative model, we compared tissue volumes, tissue distributions and accuracy at classifying non-brain intracranial tissue (arteries) and grey matter nuclei in two independent datasets (discovery N = 259, validation N = 87). Results show that compared to unimodal tissue segmentation, multimodal tissue segmentation gives more replicable volume estimations, more replicable tissue modelling, and more accurate results with regards to non-brain tissue (e.g. meninges or vessels), but only when the right model parameterization is used (i.e. 2 Gaussians per brain tissue class).

MRI

brain tissue segmentation

Gaussian Mixture Model parameterization

multimodal

SPM

Image segmentation is the generic term describing the process of partitioning an image into multiple segments, thus locating and attributing labels to different parts of an image. Many techniques exist to perform image segmentation (Gwet et al., 2018), and many have found their way into medical imaging (Dar & Padha, 2019; Lee et al., 2015). In brain imaging, tissue segmentation is a key computational element at the basis of many neuroimaging applications (González-Villà et al., 2016), for example, volume modelling and rendering, image registration and normalisation, and abnormal tissue delineation. At the heart of all those techniques is the ability to distinguish major tissue classes: grey matter, white matter, cerebrospinal fluid, meninges, and bones. The challenge in brain MRI segmentation relates to a combination of factors: natural tissue intensity variations (for instance not all grey matter voxels have the same value), partial volume effects (a given voxel can contain multiple tissue types), measurement noise (scanner readout and subjects’ motion add random variations) and spatial non-uniform intensity (i.e. because magnetic field changes over space, two voxels with the same tissue content have different values based on their location in that field). One popular approach for brain tissue segmentation is the ‘unified segmentation’ (Ashburner & Friston, 2005) approach, whereby segments are created via a Gaussian Mixture Model (GMM). The image histogram is modelled as a mixture of Gaussians representing the different tissue types and regularised using spatial priors while correcting for the bias field. While initially developed for T1 weighted images, it was extended to accommodate multiple image types, such as T1 and T2 weighted images from the same subjects (a.k.a. multispectral or multimodal segmentation). Since different MRI contrast weighting (e.g. T1w versus T2w) enhances different tissue types, multimodal segmentation is deemed better than unimodal segmentation. The choice of contrast images, depending on the underlying possible pathology, is not trivial, leading to different volume estimations and local differences (Lindig et al., 2018). For instance, in brain tumours, the combination of T1w and T2-FLAIR images might be superior to T1w and T2w images, which are more adapted for healthy individuals and maybe stroke patients. In this regard, which combination of images and how many are needed for a given pathology is still largely unknown. Here, we revisited this problem, focusing on MRI contrast image inputs (one vs two images) and the GMM parameterization because tissue classification errors can occur if the model does not adequately describe the input data. Indeed, unimodal and multimodal segmentation have been compared (Lindig et al., 2018) but without changing the underlying generative model (i.e. using the default 1 Gaussian per tissue type). Since multiple images with different contrasts are used, one can intuitively expect that more Gaussians are needed. Here, we thus quantified volumes and tissue classes systematically varying inputs and the GMM, providing recommendations on what to do and in which cases to segment the brains of individuals without noticeable abnormalities (healthy participants or psychiatric patients).

Data

The first dataset of 259 participants was prepared from 3 studies (Armand et al., 2024; Deen et al., 2018; Köhler-Forsberg et al., 2020) with data collected on two scanners (Siemens Prisma using a 32-channel head coil and Siemens Trio with Prisma fit upgrade using a 64-channel head and neck coil). All subjects were neuroanatomically normal but clinically classified as either “healthy control” or with major depressive disorder. This dataset was used as our discovery set. It had an average Male/Female ratio of 0.31 (23.94% males, 76.06% females), an average Control/Patient ratio of 1.67 (62.54% controls, 37.45% patients) and participants’ mean age was 27.29 y.o. (std 7.98). A second dataset of 87 participants from Baranger et al. (2021) was used as our validation set, which also included controls and patients suffering from depression (major depressive disorder or persistent depressive disorder, Manelis et al., 2021 - OpenNeuro ds003653). All scans were obtained from a Siemens Trio with Prisma fit upgrade using a 64-channel head coil. It had an average Male/Female sex ratio of 0.33 (25.29% males, 74.71% females), an average Control/Patient ratio of 1.23 (55.17% controls, 44.83% patients) and participants’ mean age was 28.77 y.o. (std 6.42). Further demographic information and imaging parameters for each study are presented in Table 1.

Table 1

*T1w and T2w image acquisition parameters with subjects’ demographics split per study they originate from*.
	Discovery dataset			Validation
Reference	Köhler-Forsberg et al. 2020	Deen et al., 2018	Armand et al., 2024	Manelis et al., 2021
N	141	57	61	87
Mean Age [Min Max]	26.9 [18.2 57.3]	30.1 [19.1 60.1]	25.5 [19.2 45.6]	28.77 [19.1 44.3]
Male/Female ratio	0.33	0.16	0.52	0.33
Control/Patient ratio	0.45	1	1	1.23
Scanner	Trio-Prisma-fit	Trio-Prisma-fit	Prisma	Trio-Prisma-fit
T1w parameters Resolution (mm): TR (ms): TE (ms): Flip angle (deg):	0.9x0.9x0.9 1900 2.58 9	0.9x0.9x0.9 1900 2.58 9	0.9x0.9x0.9 2000 2.58 8	0.8x0.8x0.8 2400 2.22 8
T2w parameters Resolution: TR (ms): TE (ms): Flip angle (deg):	0.9x0.9x0.9 3200 408 120	0.9x0.9x0.9 3200 408 120	0.9x0.9x0.9 3200 408 120	0.8x0.8x0.8 3200 563 120

Image Analysis

The SPM12 (version 7487, Flandin & Friston, 2008) unified segmentation algorithm (Ashburner & Friston, 2005) was used with default parameters except for the number of Gaussians set to 1 or 2 for each brain tissue class: grey matter (GM), white matter (WM), and cerebrospinal fluid (CSF). Either the T1w image only or the T1w and T2w images were used as input, thus performing the analysis four times.

Following tissue segmentation, tissue volumes were obtained for GM, WM, CSF, soft tissues and skull as part of the SPM12 pipeline. The SPM template GM, WM, CSF, and soft tissue class were summed and thresholded to keep voxels with a probability bigger than 1% of any tissue, creating a brain mask from which subjects’ tissue class data were extracted and Harell-Davis estimators of deciles computed (Harrell & Davis, 1982). Entropy and a modified version of the Dunn Index were also computed but are not reported here because entropy results did not provide more insight than other approaches, and the Dunn Index values were virtually identical across methods. All those values are available in the GitHub repository/Zenodo archive alongside the other metrics for the interested reader. To measure segmentation accuracy, tissue class probabilities were extracted for specific regions of interest: arteries and subcortical nuclei. For voxels containing arteries, the amount of GM or WM tissue should be low, while for subcortical nuclei, the amount of GM should be high, and the amount of WM and CSF should be low. The atlas of cerebral arteries (Mouches & Forkert, 2019) and subcortical nuclei (Pauli et al., 2018) were co-registered to the SPM template and the resulting image thresholded (i) to keep arteries with a diameter bigger than 0.5 millimetres and (ii) voxel nuclei with increasing probabilities (by steps of 10%).

Statistical analyses

Volume estimates

The Total Intracranial Volume (TIV) was computed by summing GM, WM, and CSF volumes. All volumes (TIV and tissue classes) were estimated using a 20% trimmed mean with 95% highest-density intervals (HDI). Note that if data are normally distributed, the trimmed mean and mean are identical, while for data with outliers, the trimmed mean is a better estimate of the population mean (Wilcox, 2012; Wilcox & Serang, 2016). We operationally defined reproducibility as an overlap of HDIs between the discovery and validation set.

Total Intracranial Volume

TIV data for the discovery set were subjected to a scanner (2 groups) x modality (uni- vs. multi-modal) x parameterisation (1 or 2 Gaussians) Hotelling T^2 repeated measures ANOVA. Significance was obtained by bootstrapping the data 1000 times under the null (i.e. after centring the data - Wilcox, 2012). For the follow-up analysis on the validation dataset, one-sample percentile t-tests on trimmed mean TIV differences were computed with a Bonferonni correction for multiple testing, testing for the effect observed in the discovery set.

Tissue volume analyses

(i) For the discovery set, scanner (2 groups) x modality (uni- vs multi-modal) x parameterisation (1 or 2 Gaussians) Hotelling T^2 repeated measures ANOVA were used for GM, WM, and CSF. Significance was obtained by bootstrapping the data 1000 times under the null hypothesis. Since interaction effects were observed, they were tested using percentile one-sample t-tests on trimmed mean differences for the validation set. All simple effects are also reported with their 95% HDI. (ii) A Gaussian Mixture Model clustering (Martinez & Martinez, 2008) of the volume differences in the 3-dimensional GM, WM, and CSF space was performed to understand how inputs and model parameters change class relations. This allowed us better to understand the effect of model parameters across subjects. Nine different models were computed using different covariance structures (clusters with (1) same variance, equal covariances, (2) same variance, unequal covariances, (3) different variance, equal covariances, (4) different variance, same volume, unequal shape, (5) different variance, unequal volume, unequal shape (6) different off-diagononals but equal covariance, (7) different off-diagononals and different orientation, (8) different off-diagononals, orientation and shape and (9) unconstrained on all aspects), and in each case, the number of clusters varied from 1 to 8. The maximum number of clusters was constrained to eight since there are only eight possible combinations for bigger/lower GM/WM/CSF. The winning models were selected using Bayesian information criteria. Post-hoc aggregation of age, sex, scanner and group belonging per cluster was also computed to check whether these hidden factors influenced the clustering algorithm.

Tissue distributions analysis

To characterise how tissue class distributions change with input and parameters, decile differences between unimodal and multimodal model inputs were analysed, a technique also known as shift function analysis (Rousselet et al., 2017), using the decile values of each tissue class (see above image analysis). Significance testing was performed using a series of one-sample t-tests with multiple comparison corrections applied among deciles for a given shift function but not across parameterisations and tissue classes.

Tissue class accuracy

(i) For voxels containing arteries, and since ‘artery’ tissue does not exist, we computed the ratio of the number of voxels below 10% and above 90% probability of a tissue class. If such a ratio is bigger than 1, it indicates that more voxels were seen as not from that tissue and, conversely, thus measuring the accuracy of the segmentation in those voxels. As previously, a scanner (2 groups) x modality (uni- vs multi-modal) x parameterisation (1 or 2 Gaussians) Hotelling T^2 repeated measures ANOVA was computed on the discovery dataset. For the follow-up analysis on the validation dataset, robust one-sample t-tests on trimmed mean differences were computed with a Bonferonni correction for multiple testing. (ii) For voxels containing GM nuclei, priors are also probabilistic, and we thus computed the GM/(WM + CSF) ratio for voxels with a-priori probabilities to be nuclei in 10%, 20%, etc .. up to 100% of subjects. For a given nucleus, the same number of voxels was used among steps (based on the smallest), but actual voxel selection was pseudo-random (constrained by the prior probability). This resulted in 9 data points per subject summarised using the trimmed mean across subjects. Curves were compared by checking the overlapping of the 99.99% confidence intervals (i.e. Bonferonni adjusted).

Volume estimates

Overall, the discovery and validation datasets gave similar volume estimates with the 95% HDIs of the trimmed means overlapping (Table 2).

Effects on TIV

The repeated measures T² ANOVA on the discovery set revealed (i) a significant effect of the modality (F(1,257) = 1748; p < .001) such that adding the T2w image decreased the TIV by ~ 37ml; (ii) a significant effect of the parameterisation (F(1,257) = 107; p < .001) such that going from 1 Gaussian to 2 Gaussians increased TIV by ~ 1 ml and (iii) no interaction (F(1,257) = 0.57; p = .41). Although there was no main scanner effect (F(1,257) = 1.61; p = .2), it interacted with segmentation parameters (scanner * modality * parameterisation, F(1,257) = 12.31; p < .001) such that adding the T2w image using 1 Gaussian or 2 Gaussians had a difference of -1.4ml for the Prisma and 0.03ml for the Prisma-fit. The main effects observed in the discovery dataset were replicated in the validation dataset, such that adding a T2w image decreased TIV by ~ 29ml, whereas adding 1 Gaussian increased TIV by ~ 3ml (p = 0.0017). Effect sizes were comparable between datasets except for a significantly larger effect of adding 1 Gaussian in the validation dataset compared to the discovery dataset (while similar in size to the discovery sub-group using the Prisma scanner only -- see confidence intervals, Table 3 and Fig. 1).

Table 2

Trimmed mean total intracranial, grey matter, white matter, and CSF volumes with 95% Highest density intervals for the discovery (N = 259) and validation (N = 87) datasets.
	Discovery set	Validation set
	TIV
T1w input, 1 Gaussian per tissue class	1442 ml [1401 1463]	1399 ml [1311 1443]
T1w input, 2 Gaussians per tissue class	1443 ml [1405 1465]	1400 ml [1327 1444]
T1w and T2w inputs, 1 Gaussian per tissue class	1404 ml [1363 1425]	1366 ml [1253 1411]
T1w and T2w inputs, 2 Gaussians per tissue class	1405 ml [1362 1427]	1373 ml [1287 1416]
	GM
T1w input, 1 Gaussian per tissue class	741 ml [722 752]	722 ml [681 744]
T1w input, 2 Gaussians per tissue class	725 ml [703 736]	708 ml [670 727]
T1w and T2w inputs, 1 Gaussian per tissue class	721 ml [701 733]	704 ml [654 726]
T1w and T2w inputs, 2 Gaussians per tissue class	731 ml [709 741]	710 ml [673 730]
	WM
T1w input, 1 Gaussian per tissue class	457 ml [435 466]	451 ml [426 466]
T1w input, 2 Gaussians per tissue class	471 ml [450 479]	464 ml [435 481]
T1w and T2w inputs, 1 Gaussian per tissue class	431 ml [415 440]	437 [415 453}
T1w and T2w inputs, 2 Gaussians per tissue class	463 ml [444 471}	462 ml [429 478]
	CSF
T1w input, 1 Gaussian per tissue class	243 ml [233 251]	218 ml [196 231]
T1w input, 2 Gaussians per tissue class	246 ml [234 253]	221 ml [196 237]
T1w and T2w inputs, 1 Gaussian per tissue class	252 ml [240 258]	218 ml [187 232]
T1w and T2w inputs, 2 Gaussians per tissue class	212 ml [201 217]	193 ml [167 205]

Table 3

TIV trimmed mean differences with 95% HDI for simple effects.
	Discovery set	Validation set
1 Gaussian per tissue class Effect of adding a T2w image	-37.3 ml [-42 -35]	-32.2 ml [-40 -28]
2 Gaussians per tissue class Effect of adding a T2w image	-36.9 ml [-40 -35]	-27.1 ml [-33 -24]
T1w as model input Effect of adding a Gaussian	+ 0.68 ml [0.3 0.89]	+ 0.68 ml [0.14 1.06]
T1w and T2w as model inputs Effect of adding a Gaussian	+ 1.19 ml [0.62 1.48]	+ 5.69 ml [3.91 6.82]

Effects on tissue types

Overall, results indicate significant modulatory effects of adding a T2w image as input and of adding a Gaussian to the model with differential effect on tissue types. Strong effect sizes were also reproduced between datasets (overlap of 95% HDI), and all changes are summarised in Table 4.

Comparable effects were observed for grey and white matter tissues. For GM, the repeated measures ANOVA revealed a main effect of input (F(1,257) = 97; p < .001) and parameterisation (F(1,257) = 23, p < .001) and an interaction (F(1,257) = 403; p < .001). Adding a T2w image reduced GM by 21ml when using a single Gaussian but increased volume estimates by 3.58ml when using 2 Gaussians. This effect was not replicated in the validation dataset (p = 0.29), for which only a reduction in volume was observed using a single Gaussian but no effect when using 2 Gaussians. For WM, the repeated measures ANOVA revealed a main effect of input (F1,257) = 1467; p < .001) and parameterisation (F(2,257) = 8970; p < .001) and an interaction (F1,257) = 1986; p < .001). Adding a T2w image reduced WM by 25ml when using a single Gaussian and 7.9ml when using 2 Gaussians. This effect was not replicated in the validation dataset (p = 0.35), for which a reduction in volume was observed using a single Gaussian but no effect when using 2 Gaussians. There were also significant scanner * modality * parameterisation interactions in the discovery dataset. For GM (F(1,257) = 21; p < .001), the effect of adding the T2w image using 1 Gaussian or 2 Gaussian has a difference of -35ml for the Prisma and − 21ml for the Prisma-fit. For WM (F(1,257) = 75; p < .001), the effect of adding a T2w image using 1 Gaussian or 2 Gaussian has a difference of -10 ml for the Prisma and − 19 ml for the Prisma-fit.

CSF showed the opposite effects as GM/WM. The repeated measure ANOVA revealed a main effect of input (F(1,257) = 102; p < .001) and parameterisation (F(2,156) = 846;p < .001) and an interaction (F(1,257) = 1687; p < .001). Adding a T2w image increased CSF by 9.2ml when using 1 Gaussian but decreased volume estimates by 32.7ml when using 2 Gaussians. This effect was not replicated in the validation dataset (p = 0.38), for which a reduction in volume was observed using 2 Gaussians but no effect when using 1 Gaussian only. In the discovery dataset, there were also significant scanner * modality (F(1,257) = 10; p < .001) and scanner * parameter effects (F(1,257) = 25; p < .001) such as adding a T2w image has a difference of -3.7ml for the Prisma and − 13.8ml for the Prisma-fit and adding a Gaussian has a difference of -24 ml for the Prisma and − 15.8ml for the Prisma-fit.

Table 4

20% GM, WM and CSF trimmed mean volume differences width 95% HDI for simple effects.
	Discovery set	Validation set
1 Gaussian per tissue class Effect of adding a T2w image	GM: -21.3 [-28.2 -18.1] WM: -25.6 [-29 -24.1] CSF: 9.2 [3.8 12.2]	GM: -19,4 [-27.9 -16.5] WM: -13,1 [-17.6 -11.4] CSF: -0.4 [-7.6 2.8]
2 Gaussians per tissue class Effect of adding a T2w image	GM: 3.5 [-1.6 5.9] WM: -7.93 [-9.3 -6.9] CSF: -32.7 [-38.9 -29.9]	GM: 0.7 [-4.3 4] WM: 0.1 [-3.4 1,3] CSF: -28.6 [-38.9 -24.2]
T1w as model input Effect of adding a Gaussian	GM: -15.6 [-17.8 -14.5 ] WM: 13.5 [12.2 14.1] CSF: 2.87 [0.32 4.1]	GM: -18.4 [-14.7 -12.9] WM: 12.5 [10.6 13.5] CSF: 2.8 [-0.9 4.5]
T1w and T2w as model inputs Effect of adding a Gaussian	GM: 8.4 [2.9 11.77] WM: 31.5 [29.6 32.5] CSF; -38.8 [-45.5 -36]	GM: 5.4 [-4.6 9.79] WM: 25.3 [22.5 26.7] CSF: -24.8 [-32 -21.8]

Importantly, the relative increases and decreases in GM/WM vs WM/CSF reported were proportional when changing parameterisation, i.e. voxels are attributed differently to GM/WM/CSF tissue classes by changing the number of Gaussian in the model. In comparison, ‘missing’ volumes observed after adding a T2w image are compensated by a re-attribution of intra-cerebral tissue to (mostly) the soft tissue class. This is illustrated in Fig. 2, where one can see how CSF using T1w images as input has more voxels all around the edges of the brain.

Tissue volume change

Focusing on change induced by adding the T2w image as input, subjects were classified by increasing/decreasing volumes, taking the three tissue types simultaneously. The classification table (Table 5) shows that the most common change is a decrease in WM volumes. When the generative model includes 1 Gaussian per tissue class, there is also a decrease in GM, while subjects are distributed among increase/decrease of CSF (~ 45% of all subjects). When the generative model includes 2 Gaussians per tissue class, WM changes are accompanied by a consistent CSF decrease. At the same time, subjects are distributed among the increase/decrease of GM (~ 33% of all subjects).

Clustering on volume differences (Fig. 3) revealed, for both the 1 Gaussian and the 2 Gaussians parameterisation, that changes in the discovery dataset can be summarised by 3 clusters with equal covariances (model 6), pulling out one small subgroup of 4 or 5 outlier subjects with substantial CSF change (-91ml and − 131ml while the CSF 95% confidence interval upper bound of the trimmed mean is + 33ml and − 7ml). The middle size cluster (N = 25 and 33, with 22 subjects in common between the 1 Gaussian and 2 Gaussian parameterisations), while lying in a different subspace than the main cluster, did not show specific univariate characteristics. Applying this model built from the discovery dataset onto the validation dataset gave satisfactory results with equivalent mean model errors (the average of subjects’ posterior probability), indicating that the model generalises. Computing new models on the validation dataset gave different results, with 8 clusters in the 1 Gaussian parameterisation and 2 clusters in the 2 Gaussian parameterisation. In this latter case, the clustering model was similar (model 7) to the discovery dataset clustering but without the ‘outlier’ cluster, showing that the segmentation did not just give similar volume estimates but also affected images across subjects in a reproducible manner. The similarity of the model 6 in the discovery dataset and model 7 in the validation dataset is striking, looking at the covariance matrix: main cluster discovery µ= [3.5–8.8 -30, diag(𝚺)= [0.22 0.028 0.3], main cluster validation (µ= [5 -2.8 -27] diag(𝚺)= [0.22 0.03 0.27]. Post-hoc class attribution did not segregate sex, group (control/patient), scanner or age in the different clusters.

Table 5

Percentages of subjects showing concurrent GM, WM and CSF volume changes caused by adding a T2w image as input for discovery (Disc) and validation (Valid) datasets when using either 1 or 2 Gaussians in the generative model.
		GM increase		GM decrease
		WM increase	WM decrease	WM increase	WM decrease
1 Gaussian	CSF increase	Disc: 0 Valid: 0	Disc: 3.47 Valid: 2.29	Disc: 0.38 Valid: 0	Disc: 66.02 Valid: 45.97
1 Gaussian	CSF decrease	Disc: 0 Valid: 0	Disc: 13.12 Valid: 3.44	Disc: 0.38 Valid: 0	Disc: 16.6 Valid: 48.27
2 Gaussians	CSF increase	Disc: 0 Valid: 0	Disc: 0 Valid: 2.29	Disc: 0 Valid: 0	Disc: 3.47 Valid: 3.44
2 Gaussians	CSF decrease	Disc: 5.40 Valid: 29.88	Disc: 56.75 Valid: 24.13	Disc: 3.47 Valid: 25.28	Disc: 30.88 Valid: 17.24

Distribution shift function analyses

Looking at how voxels were classified across the whole brain, it is essential to note that 60% of voxels were classified as having some GM, 40% with WM and 50% with CSF (i.e. many deciles are at 0). When adding the T2w image, all voxel probabilities changed across datasets, and model parameterisations and effects were reproduced between datasets (see Fig. 4). Statistical analyses showed that with 1 Gaussian per tissue class, GM and WM deciles were lower on average, i.e. adding T2w images as input leads to voxels with smaller probabilities of GM and WM, consistent with the observed lower volume estimates. With 2 Gaussians per tissue class, the lowest deciles show a lower probability of being GM or WM, but the highest deciles show a higher probability of being GM or WM. This result shows that adding T2w images as input increased the skewness of Gaussians. Because of the symmetric changes in skewness, volume estimates are virtually unchanged despite increased tissue differentiability. For CSF, adding T2w images as input consistently led to assigning lower probabilities, particularly for the last decile and when using 2 Gaussians. Again, results are consistent with volume changes. The high inter-subject variance can explain the lack of replicability of decreased CSF volumes when using a single Gaussian.

Segmentation accuracy: In analysing voxels containing major arteries, more accurate tissue segmentation corresponds to a higher-than-one ratio (i.e., a higher low than high probability of being of a given tissue). As shown in Fig. 5, ratios were above 1, indicating that the segmentation performed well. For GM, the repeated measure ANOVA showed a main effect of input (F(1,257) = 112) and parameterisation (F(1,257) = 220) and an interaction (F(1,257) = 594 - p < 0.001). Adding T2w images increases the GM ratio when using 1 Gaussian ([0.23 0.28] p = 0.001) but decreases it with 2 Gaussians ([-0.07 -0.04] p = 0.001). This interaction effect replicated (percentile bootstrap of the difference, p = 0.001), although only the 1 Gaussian case showed an increase ([0.207 0.267] vs [-0.04 0.008]) in the validation set. For WM, similar results were observed with a main effect of input (F(1,257) = 726) and parameterisation (F(1,257) = 2041) and an interaction (F(1,257) = 1341 - p = 0.001). Adding T2w images, however, always increased the WM ratio (1 Gaussian [10.3 11.7] p = 0.001; 2 Gaussians([0.74 1.33] p = 0.001). This interaction effect was replicated, although only the 1 Gaussian case showed an increase ([3.2 4.4] vs [-1.52 -0.65]) in the validation set. Finally, for CSF, a main effect of input (F(1,257) = 522) and parameterisation (F(1,257) = 500) and an interaction (F(1,257) = 1661 - p = 0.001) were observed. In this case, adding T2w images always increases the CSF ratio: 1 Gaussian [0.80 1.08] p = 0.001, 2 Gaussians [2.14 2.48] p = 0.001 and this effect was replicable (interaction p = 0.001, 1 Gaussian case [1.34 2.038], 2 Gausians [2.29 3.1]).

For completeness, it should also be reported that, in the discovery set, the scanner had no main effect on any of the tissue ratios. The scanner type nethertheless always interacted with parameters (GM input F(1,257) = 6.29 p = .01, parameterisation F(1,257) = 5.82 p = .01, input*parametrization F(1,257) = 1.2 p = .28; WM input F(1,257) = 5.45 p = .02, parameterisation F(1,257) = 1.02 p = .29, input*parametrization F(2,256) = 0.02 p = .95; CSF input F(1,257) = 1.43 p = .23, parameterisation F(1,257) = 4.34 p = .043, input*parametrization F(1,257) = 4.94 p = .034). These scanner interaction effects indicate that while observed changes stand across scanners (no main effect and replication between datasets), the amount of change observed varies differently in different tissues even when images are tightly matched, as in this dataset.

Next, we considered the classification of tissue in subcortical nuclei. First, only the nucleus accumbens (> 90% GM probability), the mammillary nucleus (~ 60%), the caudate nucleus (~ 60%), the hypothalamus (~ 58%) and the putamen (~ 54%) showed probabilities of GM above 50% on average for voxels labelled as such from the atlas. The habenular nuclei (~ 48%) and external amygdala (~ 45%) were just below average, the ventral pallidum voxels at ~ 25%, whilst other nuclei were simply not well captured (voxel GM probabilities below 10% and mostly seen as WM). Second, across all nuclei, confidence intervals of tissue ratios (GM/(WM + CSF)) mostly overlapped, allowing to conclude that results did not diverge among segmentation approaches (Fig. 6).

Post-mortem brain volume estimates (Maxeiner & Behnke, 2008) indicate TIC ranging from 1200 ml to 1600 ml, which fits the results observed for all tested segmentation approaches. As such, changing the number of inputs or model parameters does not lead to gross errors. Without a ground truth about tissues, how can we conclude that one method is better? We contend that multimodal segmentation with 2 Gaussians by tissue class for GM, WM and CSF is better in most cases because it is more replicable and more accurate to separate non ‘brain’ tissue.

While changing the number of Gaussians leads to an exchange of probability attributions between GM/WM/CSF tissues, adding T2w images reattributes probability to soft tissues mostly (~ 25ml) and even bone (~ 5ml). However, only when using 2 Gaussians per tissue class is a consistent segmentation achieved among subjects, as illustrated by clustering results, with almost all subjects showing a consistent CSF decrease. This can be explained by the increased separability of tissues revealed by the shift function analyses. Average volumes were similar using 1 or 2 Gaussians per tissue class (Tables 1 and 2) with small variations across inputs (no interaction effect for TIV and no replicable interaction effects on tissue volumes) because effects are observed mostly in the distribution tails. This result is important for voxel-based morphometric (VBM) studies, implying that (i) among published studies, reports of tissue volume estimates remain valid no matter the parameterisation chosen, (ii) within a VBM study, GM differences between groups are also valid independently of the parameterisation (the actual amount will differ, but the statistical difference would remain) but (iii) among VBM studies, correlations with behavioural and psychological variables are likely less replicable from an imaging side when using a single Gaussian per tissue class because voxel tissue attribution is not ‘stable’ from one dataset to another (leaving aside the issue of replicability and validity of behavioural measurements, see e.g. Nikolaidis et al., 2022). Such replicability issues are likely driven by many scanner-related differences (voxel sizes but also shimming, gradients, etc), as illustrated by the interaction effects in the validation dataset. On the practical side, we recommend always using 2 Gaussians, even when only using T1w images (as recommended in the SPM documentation but not set as default).

Regarding accuracy, results indicate that the SPM segmentation performs well for non-brain tissue, with small improvements in accuracy when adding T2w images. While the single Gaussian model was more accurate for all tissue types, we still recommend using 2 Gaussians per brain tissue class as the default choice because voxels containing vessels show a consistent decrease in the probability of being CSF while not decreasing GM/WM accuracy (but leading overall more consistent classification as shown above). Finally, looking at nuclei, no significant improvement was observed by adding T2w images. Good performance was observed among all segmentation approaches for the nucleus accumbens, the mammillary nucleus, the caudate nucleus, the hypothalamus, the putamen, the habenular nuclei and the external amygdala, indicating that SPM segmentation provides good GM estimates for those regions. Assuming high-resolution images are used (1mm or less), results observed across studies should therefore match irrespective of the modelling choice. For other nuclei, performance was poor, with low GM probabilities overall. In those cases, we also noted that using multimodal segmentation with 2 Gaussians was worse than in other cases. This makes sense because this model lowers the GM probabilities for the voxels in the lowest deciles and increases the WM probabilities for the higher deciles, giving ratios below 1. While those other nuclei (pallidum, globus pallidus, parabrachial, substantia nigra, red nucleus, VTA) were typically sampled with fewer voxels, this does not seem to be the reason for weak performances (for instance, habenular and mammillary nuclei also have a small number of voxel but probabilities are high). It could be that the atlas registration was not accurate enough but, more likely, imaging these structures required different sequences than ‘standard’ T2w images as demonstrated by many deep brain stimulation studies (see e.g. Isler & Bas, 2024; Nagahama et al., 2015). In general, for accurate segmentation of deep nuclei themselves (rather than tissue classification), we recommend using specialised segmentation tools (Baniasadi et al., 2022; Beliveau et al., 2021; Greve et al., 2021).

In conclusion, we showed that using multimodal segmentation is better because it increases tissue separability, leading to (i) a more stable segmentation across subjects and replicable measurements across datasets and (ii), to a lesser extent, more accurate results for non-brain intracranial tissues like vessels, but only if the right model parameterisation is used.

The authors declare no conflict of interest.

Author Contribution

CP, MC had substantial contributions to conception, design, and analysis; CP drafted the article; MC, CS and PF substantially contributed to the interpretation, revised critically for important intellectual content and approved of the version to be published.

Acknowledgments

The discovery dataset was funded by Innovation Fund Denmark (grant number: 4108-00004B), the Lundbeck Foundation (grant no R180-2014-3398 and R281-2018-131) and the Migraine Research Foundation.

Data Availability

The raw data for the discovery set are not openly available due to data privacy restrictions. However, they can be made available upon request via the Cimbi database (Knudsen et al., 2016). The raw data for the validation set are available on OpenNeuro (Markiewicz et al., 2021 - ds003653). All intermediate data (volumes and tissue class values) are shared with our analysis code on GitHub (https://github.com/CPernet/SPM_multimodal_segmentation), allowing the entire statistical analysis to be recomputed from those files. The image analysis code is also available and can be run on any valid BIDS dataset (Gorgolewski et al., 2016).

Armand, S., Langley, C., Johansen, A., Ozenne, B., Overgaard-Hansen, O., Larsen, K., Jensen, P. S., Knudsen, G. M., Sahakian, B. J., Stenbæk, D. S., & Fisher, P. M. (2024). Functional brain responses to emotional faces after three to five weeks of intake of escitalopram in healthy individuals: A double-blind, placebo-controlled randomised study. Scientific Reports, 14(1), 3149. https://doi.org/10.1038/s41598-024-51448-2
Ashburner, J., & Friston, K. J. (2005). Unified segmentation. Neuroimage, 26(3), 839–851. https://doi.org/10.1016/j.neuroimage.2005.02.018
Baniasadi, M., Petersen, M. V., Gonçalves, J., Horn, A., Vlasov, V., Hertel, F., & Husch, A. (2022). DBSegment: Fast and robust segmentation of deep brain structures considering domain generalization. Human Brain Mapping, 44(2), 762–778. https://doi.org/10.1002/hbm.26097
Baranger, D. A. A., Halchenko, Y. O., Satz, S., Ragozzino, R., Iyengar, S., Swartz, H. A., & Manelis, A. (2021). Aberrant levels of cortical myelin distinguish individuals with depressive disorders from healthy controls. NeuroImage Clinical, 32, 102790. https://doi.org/10.1016/j.nicl.2021.102790
Beliveau, V., Nørgaard, M., Birkl, C., Seppi, K., & Scherfler, C. (2021). Automated segmentation of deep brain nuclei using convolutional neural networks and susceptibility weighted imaging. Human Brain Mapping, 42(15), 4809–4822. https://doi.org/10.1002/hbm.25604
Dar, A. S., & Padha, D. (2019). Medical Image Segmentation A Review of Recent Techniques, Advancements and a Comprehensive Comparison. International Journal of Computer Sciences and Engineering, 7(7), 114–124. https://doi.org/10.26438/ijcse/v7i7.114124
Deen, M., Hansen, H. D., Hougaard, A., da Cunha-Bang, S., Nørgaard, M., Svarer, C., Keller, S. H., Thomsen, C., Ashina, M., & Knudsen, G. M. (2018). Low 5-HT1B receptor binding in the migraine brain: A PET study. Cephalalgia, 38(3), 519–527. https://doi.org/10.1177/0333102417698708
Flandin, G., & Friston, K. J. (2008). Statistical parametric mapping (SPM). Scholarpedia, 3(4), 6332. https://doi.org/10.4249/scholarpedia.6232
González-Villà, S., Oliver, A., Valverde, S., Wang, L., Zwiggelaar, R., & Lladó, X. (2016). A review on brain structures segmentation in magnetic resonance imaging. Artificial Intelligence in Medicine, 73, 45–69. https://doi.org/10.1016/j.artmed.2016.09.001
Greve, D. N., Billot, B., Cordero, D., Hoopes, A., Hoffmann, M., Dalca, A. V., Fischl, B., Iglesias, J. E., & Augustinack, J. C. (2021). A deep learning toolbox for automatic segmentation of subcortical limbic structures from MRI images. Neuroimage, 244, 118610. https://doi.org/10.1016/j.neuroimage.2021.118610
Gwet, D. L. L., Otesteanu, M., Libouga, I. O., Bitjoka, L., & Popa, G. D. (2018). A Review on Image Segmentation Techniques and Performance Measures. 12(12), 11.
Harrell, F. E., & Davis, C. E. (1982). A new distribution-free quantile estimator. Biometrika, 69(3), 635–640.
Isler, C., & Bas, G. (2024). New radiological techniques for planning of deep brain stimulation. Deep Brain Stimulation, 4, 24–28. https://doi.org/10.1016/j.jdbs.2023.12.004
Knudsen, G. M., Jensen, P. S., Erritzoe, D., Baaré, W. F. C., Ettrup, A., Fisher, P. M., Gillings, N., Hansen, H. D., Hansen, L. K., Hasselbalch, S. G., Henningsson, S., Herth, M. M., Holst, K. K., Iversen, P., Kessing, L. V., Macoveanu, J., Madsen, K. S., Mortensen, E. L., Nielsen, F. Å., & Frokjaer, V. G. (2016). The Center for Integrated Molecular Brain Imaging (Cimbi) database. Neuroimage, 124(Pt B), 1213–1219. https://doi.org/10.1016/j.neuroimage.2015.04.025
Köhler-Forsberg, K., Jorgensen, A., Dam, V. H., Stenbæk, D. S., Fisher, P. M., Ip, C. T., Ganz, M., Poulsen, H. E., Giraldi, A., Ozenne, B., Jørgensen, M. B., Knudsen, G. M., & Frokjaer, V. G. (2020). Predicting Treatment Outcome in Major Depressive Disorder Using Serotonin 4 Receptor PET Brain Imaging, Functional MRI, Cognitive-, EEG-Based, and Peripheral Biomarkers: A NeuroPharm Open Label Clinical Trial Protocol. Frontiers in Psychiatry, 11, 641. https://doi.org/10.3389/fpsyt.2020.00641
Lee, L. K., Liew, S. C., & Thong, W. J. (2015). A Review of Image Segmentation Methodologies in Medical Image. In H. A. Sulaiman, M. A. Othman, M. F. I. Othman, Y. A. Rahim, & N. C. Pee (Eds.), Advanced Computer and Communication Engineering Technology (Vol. 315, pp. 1069–1080). Springer International Publishing. https://doi.org/10.1007/978-3-319-07674-4_99
Lindig, T., Kotikalapudi, R., Schweikardt, D., Martin, P., Bender, F., Klose, U., Ernemann, U., Focke, N. K., & Bender, B. (2018). Evaluation of multimodal segmentation based on 3D T1-, T2- and FLAIR-weighted images—The difficulty of choosing. Neuroimage, 170, 210–221. https://doi.org/10.1016/j.neuroimage.2017.02.016
Manelis, A., Baranger, D. A. A., Satz, S., Ragozzino, R., Lyengar, S., Swartz, H. A., & Halchenko, Y. O. (2021). Cortical myelin measured by the T1w/T2w ratio in individuals with depressive disorders and healthy controls. [BIDS]. OpenNeuro. https://doi.org/10.18112/openneuro.ds003653.v1.0.0
Markiewicz, C. J., Gorgolewski, K. J., Feingold, F., Blair, R., Halchenko, Y. O., Miller, E., Hardcastle, N., Wexler, J., Esteban, O., Goncavles, M., Jwa, A., & Poldrack, R. (2021). The OpenNeuro resource for sharing of neuroscience data. eLife, 10, e71774. https://doi.org/10.7554/eLife.71774
Maxeiner, H., & Behnke, M. (2008). Intracranial volume, brain volume, reserve volume and morphological signs of increased intracranial pressure – A post-mortem analysis. Legal Medicine, 10(6), 293–300. https://doi.org/10.1016/j.legalmed.2008.04.001
Mouches, P., & Forkert, N. D. (2019). A statistical atlas of cerebral arteries generated using multi-center MRA datasets from healthy subjects. Scientific Data, 6(1). Article 1. https://doi.org/10.1038/s41597-019-0034-5
Nagahama, H., Suzuki, K., Shonai, T., Aratani, K., Sakurai, Y., Nakamura, M., & Sakata, M. (2015). Comparison of magnetic resonance imaging sequences for depicting the subthalamic nucleus for deep brain stimulation. Radiological Physics and Technology, 8(1), 30–35. https://doi.org/10.1007/s12194-014-0283-0
Nikolaidis, A., Chen, A. A., He, X., Shinohara, R., Vogelstein, J., Milham, M., & Shou, H. (2022). Suboptimal phenotypic reliability impedes reproducible human neuroscience (p. 2022.07.22.501193). bioRxiv. https://doi.org/10.1101/2022.07.22.501193
Pauli, W. M., Nili, A. N., & Tyszka, J. M. (2018). A high-resolution probabilistic in vivo atlas of human subcortical brain nuclei. Scientific Data, 5(1), 180063. https://doi.org/10.1038/sdata.2018.63
Rousselet, G. A., Pernet, C. R., & Wilcox, R. R. (2017). Beyond differences in means: Robust graphical methods to compare two groups in neuroscience. European Journal of Neuroscience, 1738–1748. https://doi.org/10.1111/ejn.13610
Wilcox, R. R. (2012). Introduction to robust estimation and hypothesis testing (3rd ed.). Academic.
Wilcox, R. R., & Serang, S. (2016). Hypothesis Testing, p Values, Confidence Intervals, Measures of Effect Size, and Bayesian Methods in Light of Modern Robust Techniques. Educational and Psychological Measurement. https://doi.org/10.1177/0013164416667983

No competing interests reported.

Download PDF

Editorial decision: Revision requested
21 Aug, 2024
Reviews received at journal
21 Aug, 2024
Reviews received at journal
18 Aug, 2024
Reviewers agreed at journal
09 Aug, 2024
Reviewers agreed at journal
08 Aug, 2024
Reviewers invited by journal
08 Aug, 2024
Editor assigned by journal
08 Aug, 2024
Submission checks completed at journal
08 Aug, 2024
First submitted to journal
06 Aug, 2024

You are reading this latest preprint version

Multimodal Tissue Segmentation is better

Status:

Version 1

Abstract

Figures

Introduction

Method

Data

Image Analysis

Statistical analyses

Results

Discussion

Declarations

Author Contribution

Acknowledgments

Data Availability

References

Additional Declarations

Status:

Version 1