Clustering of fMRI data: the elusive optimal number of clusters

Fuzzy clustering

Our clustering method was based on the popular fuzzy c-mean (FCM) algorithm (Bezdek, 1981; Bezdek et al., 1997). In the context of fMRI, the FCM algorithm can segregate or cluster n brain voxels (feature vectors) into c expected clusters (c ≥ 2). Each voxel i is a vector X_i of p properties (e.g., number of collected volumes or scans). Each cluster j is characterised by a centroid V_j, that represents its characteristic timecourse (prototype). The resemblance between each voxel i and each centroid V_j is assessed by the distance D_ij between X_i and V_j. The degree of membership U_ij is calculated for each voxel i by comparing D_ij for each cluster j to all other clusters.

In brief, the standard FCM algorithm iteratively minimises the following objective function J_m: (1)${\displaystyle J_m=\sum _{i=1}^n\sum _{j=1}^c{U}_{ij}^m\cdot D_{ij}^2}$ where “m” is the degree of fuzziness.

Degrees of membership U and centroids V are updated as following:

(2)${\displaystyle U_{ij}=\frac{1}{\sum _{k=1}^c{\left(\frac{D_{ij}}{D_{ik}}\right)}^{2∕m-1}}}$ (3)${\displaystyle V_j=\frac{\sum _{i=1}^n{U}_{ij}^m\cdot X_i}{\sum _{i=1}^n{U}_{ij}^m}.}$

Optimal clustering depends on the choice of the similarity D, the degree of fuzziness m and the optimal number of clusters c_opt, as detailed below.

Similarity measure D

Here I used a modified version of the hyperbolic correlation distance proposed previously by Golay et al. (1998). In their work, D was defined as (Golay et al., 1998): (4)${\displaystyle D_{ij}=\frac{1-C{C}_{ij}}{1+C{C}_{ij}}.}$

Where CC_ij is the Pearson correlation coefficient between X_i and V_j.

Here, a modified version of D was used: (5)${\displaystyle D_{ij}=\frac{\sqrt{\left|C{C}_{ij}\right|}-C{C}_{ij}}{\sqrt{\left|C{C}_{ij}\right|}+C{C}_{ij}}.}$

This new formula uses the square root function, a monotonically increasing function over x > 0 that satisfies the following inequality: $\sqrt{x}\ge x$, for x ∈ [0, 1]. The rationale here was to increase the difference (i.e., discrimination power) between relatively close correlation values in particular between mid and high correlations (cf. Fig. S1).

Optimal number of clusters

A good and robust clustering should yield compact and well-separated clusters. This is assumed to be the case when the number of clusters reaches an optimal value c_opt. The exact c_opt value is however unknown in fMRI data. Previous reports have suggested that c_opt can be found within the interval [2, $\sqrt{n}$] (Zahid, Limouri & Essaid, 1999); however the exact c_opt can only be estimated empirically. Typically, FCM is repeated several times with different c values (i.e., equivalent to an unsupervised fuzzy clustering analysis Fadili et al., 2001) and the c value that optimises a given criterion, here a given CV index, is considered as the optimal c_opt, and that criterion is typically defined as a trade-off between compactness and separation.

Before introducing the different CV indices used here, it might be helpful to define the core measures of compactness and separation using unified mathematical notations. These measures can be seen as building blocks that can be combined into different CV indices. Ultimately, the definition of those measures would help appreciate the inherent links (or similarity) between previously suggested CV indices, before introducing the rationale of the new CV index.

Compactness and separation measures

Two core quantities, noted n_m,j and σ_m,j, were defined as following:

(6)${\displaystyle n_{m,j}=\sum _{i=1}^n{U}_{ij}^m}$ (7)${\displaystyle {\sigma }_{m,j}=\sum _{i=1}^n{U}_{ij}^m\cdot D_{{}_{ij}}^2.}$

The measures n_1,j and n_2,j represent the fuzzy cardinality and the fuzzy partition of cluster j respectively. The quantity σ_m,j denotes the fuzzy variation of cluster j, though other studies have instead used σ_1,j as a measure of fuzzy variation (e.g., Gath & Geva, 1989; Rezaee, Lelieveldt & Reider, 1998; Sun, Wang & Jiang, 2004).

Those core quantities can then be combined into different forms to give away different measures of fuzzy compactness (cohesiveness) for a given c-partition. Using similar notation as previous studies, quantities called π_m,1 (Bensaid et al., 1996; Zahid et al., 1999), π_m,m (Bouguessa, Wang & Sun, 2006), and FC (Fadili et al., 2001; Zahid et al., 1999) were computed as following:

(8)${\displaystyle {\pi }_{m,1}=\sum _{j=1}^c\frac{{\sigma }_{m,j}}{n_{1,j}}}$ (9)${\displaystyle {\pi }_{m,m}=\sum _{j=1}^c\frac{{\sigma }_{m,j}}{n_{m,j}}}$ (10)${\displaystyle FC=\frac{\sum _{i=1}^n{\left(\underset{j}{max}\left(U_{ij}\right)\right)}^2}{\sum _{i=1}^n\underset{j}{max}\left(U_{ij}\right)}.}$

Likewise, the fuzzy separation (isolation) between clusters was previously estimated with several fuzzy separations quantities called K_m (Fukuyama & Sugeno, 1989), FS (Fadili et al., 2001; Zahid et al., 1999), S (Zahid et al., 1999) and SS (Rezaee, Lelieveldt & Reider, 1998):

(11)${\displaystyle K_m=\sum _{j=1}^c{n}_{m,j}\cdot {∥V_j-\overline{X}∥}^2}$ (12)${\displaystyle FS=\sum _{j=1}^{c-1}\sum _{k=1}^{c-j}\frac{\sum _{i=1}^n{\left(min\left(U_{ij},U_{i,k+j}\right)\right)}^2}{\sum _{i=1}^{n}min\left(U_{ij},U_{i,k+j}\right)}}$ (13)${\displaystyle S=\frac{1}{c}\sum _{j=1}^c{∥V_j-\overline{X}∥}^2}$ (14)${\displaystyle SS=\sum _{j=1}^c\frac{1}{\sum _{k=1}^c∥V_j-V_k∥}}$ where $\overline{X}$ stands for the global mean of the whole data.

Interestingly, the ratio FS/FC (i.e., separation divided by compactness) is known as the fuzzy overlap (FO) coefficient (see Fadili et al., 2000 for more details).

Furthermore, different measures of between-centroid distance have been proposed, including the minimum distance V_dmin (e.g., Schwämmle & Jensen, 2010; Xie & Beni, 1991), the maximum distance V_dmax (e.g., Rezaee, Lelieveldt & Reider, 1998), and the minimum distance V_dmin,j between a cluster j and the remaining clusters (Wu & Yang, 2005):

(15)${\displaystyle V_{dmin}=\underset{j,k}{min}\left(∥V_j-V_k∥\right)}$ (16)${\displaystyle V_{dmax}=\underset{j,k}{max}\left(∥V_j-V_k∥\right)}$ (17)${\displaystyle V_{dmin,j}=\underset{k\ne j}{min}\left(∥V_j-V_k∥\right).}$

These measures, based on the distance between estimated centroids, can be seen as alternative separation measures. They can be handy when the clustering is showing redundant clusters.

This section introduces two new measures of separation and discrimination between voxels by combining different measures of fuzzy cardinality and variation (cf. Eqs. (6) and (7)): a fuzzy intra-cluster (ID_intra) dissimilarity coefficient and an inter-cluster (ID_inter) dissimilarity coefficient:

(18)${\displaystyle I{D}_{\mathrm{int}ra}=\underset{j}{max}\left(\frac{n-n_{1,j}}{n_{1,j}}\cdot \frac{{\sigma }_{1,j}}{\sum _{k=1,k\ne j}^c{\sigma }_{1,k}}\right)}$ (19)${\displaystyle I{D}_{\mathrm{int}er}=\underset{j}{min}\left(\frac{\underset{k,k\ne j}{min}\left({\sigma }_{1,k}\right)}{{\sigma }_{1,j}}\right).}$

Small ID_intra values would indicate that, across all clusters, voxels that are close to a given cluster are well-isolated from voxels that are far from that cluster, whereas high ID_inter values indicate well-discriminated voxels (i.e., small fuzzy overlap between clusters). Our initial tests with noisy simulated data showed the need to define new separation measures that are robust to noise and can handle c-partitions with higher number of clusters, hence the new definitions in Eqs. (18) and (19).

Cluster validity measures

There are many CV indices in the literature (probably more than 50 indices), hence it is beyond the scope of this study to test all of them. In a preliminary analysis (results not shown here), about 20 selected CV indices were first tested on several simulated datasets (as defined in Bezdek & Pal, 1998; Bouguessa, Wang & Sun, 2006; Dave, 1996; Fukuyama & Sugeno, 1989; Geva et al., 2000; Kim, Park & Park, 2001; Kim, Lee & Lee, 2003; Kim & Ramakrishna, 2005; Kwon, 1998; Pakhira, Bandyopadhyay & Maulik, 2004; Pakhira, Bandyopadhyay & Maulik, 2005; Pal & Bezdek, 1995; Rezaee, Lelieveldt & Reider, 1998; Rhee & Oh, 1996; Sun, Wang & Jiang, 2004; Tsekouras & Sarimveis, 2004; Wu & Yang, 2005; Xie & Beni, 1991; Yu & Li, 2006; Zahid et al., 1999; Zahid, Limouri & Essaid, 1999). These CV indices were selected from earlier studies (for a similar rationale, see recent comparison study Zhou et al., 2014), and many of them are well-established indices. Some of these CV indices have been used in previous fMRI studies. More recent CV indices (e.g., see He, Tan & Fujimoto, 2016; Hu et al., 2011; Lin et al., 2016; Ren et al., 2016; Rezaee, 2010; Yang et al., 2018; Zhang et al., 2014) were not explicitly tested here.

From this preliminary analysis, seven CV indices (out of twenty) were selected according to the following four criteria: CV indices should (i) combine both measures of separation and compactness; (ii) not suffer from monotonic dependency with the number of expected clusters; (iii) not necessitate the categorisation or the binarisation of U (i.e., crisp degrees of membership) during CV computation; (iv) be fast to compute when n is expected to be very high (e.g., hundreds of thousands of voxels in the context of fMRI data). The seven selected CV indices that satisfied the different criteria are described below and listed in Table 1.

(1)- The Rezaee-Lelieveldt-Reider index CV_RLR (Rezaee, Lelieveldt & Reider, 1998): (20)${\displaystyle C{V}_{RLR}=\frac{\sum _{j=1}^c{\sigma }_{1,j}}{c\cdot ∥{\sigma }_X∥}+\frac{1}{\alpha }\cdot \left(\frac{V_{dmax}\cdot SS}{V_{dmin}}\right).}$

The constant α is a weighting constant and σ_X is the variance of the whole data set. The best c-partition is obtained by minimising CV_RLR with respect to the number of clusters c. In the original definition of CV_RLR, the constant α was set to 1; however, here α was set to the value of $\frac{V_{dmax}}{V_{dmin}}\cdot SS$ at the maximum number of clusters (c_max) as suggested previously (Sun, Wang & Jiang, 2004).

(2)- The Zahid-Limouri-Essaid index CV_ZLE (Zahid et al., 1999; Zahid, Limouri & Essaid, 1999): (21)${\displaystyle C{V}_{ZLE}=\alpha \cdot \left(\frac{S}{{\pi }_{m,1}}\right)-\frac{FS}{FC}.}$

The constant α is independent from c and was introduced here as a scaling factor to take into account the difference in values between the two subtracted quantities. The constant α was set here to the value of the fuzzy overlap (FS/ FC) at c = c_max (note that in the original paper of Zahid et al., α was set equal to 1). The best c-partition is obtained by maximising CV_ZLE with respect to c. This CV_ZLE index has previously been used for fMRI analysis (Fadili et al., 2001).

Note that the ratio $\left(\frac{S}{{\pi }_{m,1}}\right)$ in Eq. (21) is also known as the Pal-Bezdek cluster validity index (Pal & Bezdek, 1995).

(3)- Among several CV indices suggested by Geva and colleagues (Geva et al., 2000), the invariant index CV_GV was selected here to measure the ratio of the between-cluster scatter matrix to the within-cluster scatter matrix (Geva et al., 2000): (22)${\displaystyle C{V}_{GV}=\frac{K_1}{c^2\cdot J_1}.}$

The normalisation with the number of clusters c minimise the monotonically increase of CV_GV when c increased. This index should be maximal at the optimal c-partition.

(4)- The Kim-Park index, noted CV_KP (Kim, Park & Park, 2001): (23)${\displaystyle C{V}_{KP}=\frac{{\pi }_{1,1}}{c}+\frac{1}{\alpha }\cdot \left(\frac{c}{V_{dmin}}\right).}$

The best c-partition is obtained by minimising the index CV_KP with respect to the number of clusters c. This index has previously been used for fMRI analysis (Moller et al., 2002).

(5)- The Pakhira-Bandyopadhyay-Maulik index CV_PBM (Pakhira, Bandyopadhyay & Maulik, 2004; Pakhira, Bandyopadhyay & Maulik, 2005): (24)${\displaystyle C{V}_{PBM}=\frac{\alpha }{c}\cdot \frac{V_{dmax}}{J_m}.}$

With α as a constant term (e.g., α was set here to n). The best c-partition is obtained by maximising CV_PBM with respect to the number of clusters c.

(6)- The Wu-Yang index CV_WY (Wu & Yang, 2005): (25)${\displaystyle C{V}_{WY}=\sum _{j=1}^c\left(\frac{n_{1,j}}{\underset{j}{max}\left(n_{1,j}\right)}-exp\left(-\frac{V_{{}_{dmin,j}}^2}{S}\right)\right).}$

This index compared the fuzzy partition of each cluster to its exponential separation, with −c < CV_WY < c, and CV_WY is maximal at c_opt.

(7)- The Bouguessa-Wang-Sun index CV_BWS index (Bouguessa, Wang & Sun, 2006): (26)${\displaystyle C{V}_{BWS}=\frac{K_m}{{\pi }_{m,m}}.}$

This index CV_BWS should be maximised with respect to c.

(8)- Our new CV index, noted CV_new, combined different measures of compactness and separation as following: (27)${\displaystyle C{V}_{new}=K_m\cdot \left(\frac{I{D}_{\mathrm{int}er}}{I{D}_{\mathrm{int}ra}}\right)\cdot \left(\frac{FC}{J_1}\right).}$

The best c-partition should maximise CV_new. The rationale behind incorporating those specific compactness and separation measures (ID_inter, K_m, FC, ID_intra, J₁) in the definition of CV_new is illustrated below with simulated (noisy) datasets.

Simulated data

Twenty-two simulated datasets were generated as following. First, a fixed number c of time-courses with p datapoints (p = 100) were generated from a unit normal distribution (mean = 0, σ = 1). The Pearson correlation between these c time-courses was less than 0.1 for all simulated datasets. Second, each time-course was replicated r_j times, with j = 1…c_opt, and ${\sum }_{j=1}^c{r}_j=n$ (where n is the total number of voxels, set here to 1,000). Third, n random timecourses with p datapoints, generated from a normal distribution (mean = 0) but with variable noise levels (σ = 1 or 4) were added to the replicated time-courses. This would help to test the robustness of FCM at different noise levels (for a similar rational see Kim, Park & Park, 2001; Wang & Zhang, 2007) and to monitor the behaviour of the different CV indices when the fuzzy compactness of clusters became very low (i.e., high intra-class dissimilarity in noisy data). This procedure generated a dataset X of n voxels, each with p datapoints, with known and fixed numbers of classes. Multidimensional scaling (MDS) tools were used to visualise the simulated c clusters.

Specifically, the following 22 datasets were generated: (i) a single-cluster dataset (noted 1-cluster; e.g., the ‘null’ case, see Tibshirani, Walther & Hastie, 2001) with highly similar voxels (c_opt = 1; Fig. S2A); (ii) a dataset without any obvious structure (noted n-cluster data; c_opt near to n; Fig. S2B), see Suleman (2017); (iii) ten datasets with known number of clusters c_opt varying from 2 to 11 and low noise level (σ = 1, see illustration in Fig. S2C with c_opt = 3); (iv) ten datasets with a known number of clusters c_opt varying from 2 to 11 and high noise level (σ = 4, see illustration in Fig. S2D with c_opt = 3).

All simulated datasets were clustered by FCM with c (i.e., number of expected clusters) varying between c_min = 2 and c_max = 19. All analyses were carried out with homemade Matlab-based scripts (MathWorks, Natick, MA, USA).

Real fMRI data

Real data consisted of single subject fMRI data with a block paradigm design (freely available at: http://www.fil.ion.ucl.ac.uk/spm/data/auditory.html). The block paradigm consisted of alternated epochs between rest and auditory stimulation. 96 volumes were acquired on a modified 2T Siemens MAGNETOM Vision system (TR = 7s, 64 contiguous slices). To avoid T1 effects in the initial scans, the first 12 scans were discarded, leaving 84 scans for further analysis (p = 84). The data were realigned, normalised (voxel size 2 × 2 × 2) and smoothed (FWHM = 6 × 6 × 6 mm). This dataset was selected because it has been used in many previous studies with clustering techniques including FCM (e.g., Gu et al., 2005; Lu, Jiang & Zang, 2004). FCM was applied on this real fMRI dataset with c varying between c_min = 2 and c_max = 39. To identify relevant FCM cluster(s) with activated auditory regions, the centroids (prototype) V_j (j = 1…c) were correlated with the experimental block design (Bandettini et al., 1993).

To appreciate the distribution of brain regions’ sizes in each c-partition, a morphological granulometry was applied to all identified clusters after binarization (Soille, 2003). This analysis estimated the size of each spatially distinct region or blob (26-connected neighbourhood) for a given crisp FCM partition. Given that each voxel belongs to all clusters at different degrees of membership (cf. U_ij in Eq. (2)), the threshold was set to 0.5 so that each voxel belongs maximally to one cluster. Practically, for a given c-partition and for each binary cluster j (j = 1…c), the size of each region as well as the number of isolated voxels (i.e., single-voxel regions) were calculated.

This dataset was also analysed with SPM12 software package (Wellcome Trust Centre for Neuroimaging, London UK; http://www.fil.ion.ucl.ac.uk/spm/) using standard procedures. This allowed auditory activations to be identified using model-based methods.

Degree of fuzziness m

The degree of fuzziness m might influence the output of clustering (e.g., Bezdek, 1981; Fadili et al., 2000; Fadili et al., 2001; Krishnapuram & Keller, 1993; Selim & Ismail, 1986; Yu, Cheng & Huang, 2004): when m tends to 1 the classification becomes crisp and U_ij takes the value 0 (voxel i is not a member of cluster j) or 1 (voxel i belongs to cluster j) but when m tends to +∞ the classification is purely fuzzy (U_ij is near to 1/c). The optimal value of m may depend on the characteristics of the data. Previous empirical work approximated m by a nonlinear function of the dimensions of the data (n and p); for example, Eq. (5) of Schwämmle & Jensen (2010, Page 2845) yields m values of 1.044 and 1.019 for our artificial and real datasets respectively. However, these estimated values are too low compared to typical m values encountered in neuroimaging studies. Previous studies have explored the influence of m on the computation of CV indices (e.g., Zhou, Fu & Yang, 2014), and they found better clustering results with m between 1.2 and 2.5 for fMRI data (Fadili et al., 2000; Fadili et al., 2001; Moller et al., 2002; Smolders et al., 2007).

More specifically, there are two issues to be considered when selecting m during the computation of CV indices. First, several CV indices became inadequate with hard c-partitions (i.e., m tends to 1). Specifically, any measures that are based exclusively on the distribution of U values (e.g., FC, FS) would artificially reach their optimal values independently from the number of clusters c. Second, according to Eq. (3), centroids become close to the mean of the whole data set $\overline{X}$ when m tends towards +∞. In other words, the c clusters would have comparable fuzzy cardinality values (i.e., Eq. (6)) for larger m values, which may be problematic when some clusters are expected to contain a small number of voxels (see illustration in Fig. S3); for more details see (Selim & Ismail, 1986; Tsekouras & Sarimveis, 2004; Yu, Cheng & Huang, 2004). This issue is particularly critical when analysing task-related fMRI data because activated voxels are expected to represent a small fraction of the whole brain.

Here, m was held to 1.5 throughout this study.

Voxel selection and the ill-balanced dataset problem in fMRI

One important issue during the clustering of fMRI datasets is the selection of the relevant n voxels. Because the number of activated voxels is small (i.e., a few percent) compared to the total number of voxels in a typical whole-brain fMRI dataset, previous studies have suggested different approaches to overcome this ‘ill-balanced’ data problem. For instance, FCM can be limited to relevant voxels within the gray matter, in specific anatomical brain regions, or to voxels with some kind of task-related effects (e.g., see Fadili et al., 2000; Goutte et al., 1999; Gu et al., 2005; Lee et al., 2012; Moller et al., 2002; Seghier, Friston & Price, 2007). Voxel selection might be useful for: (i) reducing the high dimensionality of the problem and improving both computational robustness and speed; (ii) minimising the influence of redundant voxels; and (iii) increasing the accuracy of the clustering by focusing mainly on meaningful voxels. However, the author preferred here to include all brain voxels so that the robustness of the different CV indices can be appreciated when noisy voxels (voxels with no effect of interest) and artefacts are present. FCM was thus applied to all voxels of the real fMRI dataset, yielding a total number of voxels n = 227,716.

FCM convergence and the initialisation problem

Depending on the initialisation of the degrees of membership U (Bezdek, 1981), the FCM algorithm may converge to different c-partitions (e.g., local minima). This problem of initialisation may lead to spurious c-partitions (Moller et al., 2002) when using CV indices. One possible solution is to repeat the FCM algorithm on the same dataset with several different random initialisations (e.g., Moller et al., 2002; Pena, Lozano & Larranaga, 1999), with the expectation that it is unlikely that different starting conditions will lead to the same local minima. Accordingly, for each c value, the FCM algorithm was re-run on the real fMRI dataset ten times with random initialisations (for a similar procedure see Chuang et al., 1999).

FCM on simulated data

CV indices on data with known numbers of clusters

The different measures of compactness and separation are shown in Fig. 2A for the 7-clusters data set. Several measures showed different values over the number of clusters as compared to the clustering of the 1-cluster and n-cluster datasets. For instance, the coefficient J₁ decreased in the interval c = 2 to c = 7, consistent with the fact that data can be clustered further (as seen for the n-cluster data); then it reached a plateau for higher number of classes, consistent with the fact that the data cannot be segregated any further (as the case of the 1-cluster data). The limit between the two behaviours was indeed at the true number of clusters (c = 7). This observation is valid for the other measures of compactness (e.g., FC, π_m,1, ID_intra) and separation (e.g., FS, ID_inter, S, K_m). When the data became noisy, some measures were less sensitive to the structure of the data (i.e., the presence of seven clusters). As illustrated in Fig. 2B, fuzzy separation FS and compactness π_m,1 showed comparable behaviour as in the clustering of the n-cluster dataset, which reflects the influence of noisy distant points (low within-cluster compactness and between-cluster separation). Interestingly, in addition to V_dmin, quantities ID_intra, K_m and J₁ were more robust to noise and showed high discriminability with an optimal value around the expected number of classes (Fig. 2B). This observation further motivated their inclusion in the definition of the new CV index.

Figure 3 illustrates all CV indices for the 3-cluster, 7-cluster, and 11-cluster datasets with low noise level (σ = 1). All CV indices indicated the best c-partition for the expected number of clusters (maximum value for CV_ZLE, CV_GV, CV_PBM, CV_WY, CV_BWS, CV_new; minimum value for CV_RLR and CV_KP). Note that the new index CV_new is highly discriminative in pointing to the optimal c-partition. When the data became noisy (σ = 4), all CV indices, except CV_BWS and CV_new, failed to indicate the optimal c-partition (Fig. 4). However, for data with higher c_opt (e.g., c_opt > 9), only the new index CV_new identified the true number of clusters, albeit with lower discriminability (e.g., compare Figs. 3B to 4B).

An ad hoc analysis was conducted to monitor the behaviour of CV_new over different degrees of fuzziness m (m varying between 1.2 and 2.5), for a similar rationale see (Schwämmle & Jensen, 2010). This analysis showed that CV_new correctly identified the true number of clusters c_opt in almost all simulated datasets for m ∈ [1.2, 2.5], except for datasets with both high noise level (σ = 4) and high number of true clusters (c_opt > 9) where CV_new failed to identify c_opt when m ≥ 2 (i.e., CV_new underestimated c_opt at higher m values including the popular value of m =2). This ad hoc analysis confirmed the initial choice of m = 1.5.

FCM real fMRI data

As expected, the number of iterations for the convergence of the FCM algorithm varied across the 10 different initialisations. However, for a given c value and across the ten runs, the obtained c-partitions were very similar and the function J_m (Eq. (1)) reached the same minimum value (except for c values between 12 and 15 where one initialisation reached a different minimal J_m value compared to the other nine initialisations).

Identified clusters

Figure 5 plots the different coefficients and CV indices against the number of expected clusters c varying from 2 to 39. Measures such as ID_inter, and V_dmin showed an interesting pattern when c increased, with high and decreasing values for small number of clusters (c < 10) and low and fixed values when c increased (a comparable behaviour was also seen for FC). This mirrored their behaviour during the clustering of the 1-cluster and n-cluster datasets. The change in the nature of the dependency occurred around c = 13, indicating the maximum c value that ensured different centroids V. For a number of expected clusters bigger than 13, the c-partition contained a few redundant classes (identical centroids V). However, for c < 13 clusters, although the obtained classes were compact (e.g., high FC values), the separation between clusters was not optimal (see for instance V_dmax, S, and K_m). More specifically, the fuzzy separation measures S and K_m showed optimal values for higher numbers of expected clusters at c larger than 17 clusters. At this range, the c-partition contained at least three similar centroids.

Figure 5B illustrates the dependency of different CV indices with c. Some CV indices (e.g., CV_RLR and CV_KP) showed optimal values for low c values (maximal fuzzy compactness), whereas other CV indices (e.g., CV_ZLE, CV_BWS, CV_GV, and CV_PBM) showed optimal values at an intermediate number of expected clusters (i.e., maximal fuzzy separation). Interestingly, the new index CV_new went through different phases (i.e., different plateaus), depending on the weight of fuzzy separation and compactness (a change of behaviour visible at c = 15). The new index CV_new reached its maximum value at c = 24 clusters, ensuring a good compromise between separation and compactness of the c-partition of this real fMRI dataset.

The results of the morphological granulometry at U > 0.5 are illustrated in Fig. 6. As expected, the size of very large regions tend to decrease with the number of expected clusters, as large regions were subdivided further into smaller regions at higher c values. Interestingly, for each c-partition, the total number of single-voxel regions over all clusters was less than 0.04% of n (Fig. 6). Given the spatial smoothness of the fMRI data, there was no cluster containing exclusively single-voxel regions.

Figure 7 illustrates all obtained clusters for c-partitions with low fuzzy separation (c = 8), without redundant clusters (c = 13), at high fuzzy separation (c = 18), and at the optimal c value that maximised CV_new (c = 24). Identified voxels within the auditory cortex (i.e., voxels of interest) are shown in the first axial slice of each c-partition. Voxels in the auditory cortex were grouped with those in the occipital lobe at small c values (c = 8), but they became clearly segregated at larger c values (e.g., c = 18 and c = 24). Interestingly, identified voxels within the auditory cortex in the c-partition with 24 clusters were remarkably similar to those identified with model-based SPM methods (e.g., SPM map at p < 0.05 FWE-corrected, Fig. 8). Last but not least, the centroid of the relevant cluster with activations in auditory regions (Cluster “1” of the 24-partitions in Fig. 7) was strongly correlated with the experimental block design (r = 0.7, p < 0.001).