HD issues in clustering

Cluster analysis is one of the main data analysis method. It aims at estimating a K groups partition z = (z 1 , . . . , zn) of a data set x = (x 1 , . . . , xn) 2 X d Model-based clustering allows to reformulate cluster analysis as a well-posed estimation problem both for the partition z and for the number K of groups [McLachlan and Basford,88] Model-based clustering has led practical successes, even in the challenging HD-d case [Biernacki and Maugis,17] But less works on model-based clustering dedicated to the HD-m j case

Have now a look at

⇡ 1 = ⇡ 2 = 1 2 , X 1 |z 11 = 1 ⇠ N d (0, ⌃), X 1 |z 12 = 1 ⇠ N d (1, ⌃)
Each variable provides equal and own separation information HD-m j : limit of the classical Latent Class Model (LCM)

When X corresponds to d categorical variables (frequent situation), the jth one having m j response levels,the standard latent class model (LCM, [Goodman,74]) is used

p(x i ; ✓) = K X k=1 ⇡ k d Y j=1 m j Y h=1 (↵ jh k ) x jh i , where ✓ = (⇡, ↵), ⇡ = (⇡ k ), ↵ = (↵ k )
When some m j s are large another kind of dimensionality curse appears

The standard latent class model is a↵ected by such a curse of dimension paradigm since the number of free parameters if equal to K 1 + K P d j=1 (m j 1) An early work for parsimony in HD-m j [Lebret et al.,15] proposes four parsimonious versions of LCM. They correspond to an extension of the parameterization of Bernoulli distributions used by [Celeux and Govaert,91] for clustering and also by [Aitchinson and Aitken,76] for kernel discriminant analysis.

The basic idea is

↵ j k = (↵ j1 k , . . . , ↵ jm j k ) = ( j k , . . . , j k , j k , j k , . . . , j k ) with j k > j k . Since Pm j h=1 ↵ jh k = 1, then (m j 1) j k + j k = 1 and, thus, j k = (1 j k )/(m j 1).
It leads to a very parsimonious (only K 1 + Kd parameters) but poor flexible model. We propose to extend this idea by a level partitioning (here an example with 3 level clusters)

↵ j k = (↵ j1 k , ↵ j2 k , ↵ j3 k , ↵ j4 k , ↵ j5 k ) = ( j 1 k , j 2 k , j 1 k , j 3 k , j 2

Related EM algorithm

As a consequence of the preeceeding remark, estimating # can be performed by any classical EM algorithm [Dempster et al.,77] for LCM implemented in any classical existing package (for instance Rmixmod [Lebret et al.,15] or RMixtComp1 ).

Starting from # 0 = (⇡ 0 , 0 ), iteration r > 0 of this EM is expressed as:

E step: Conditional probability that individual i arose from cluster k t ik (# r ) = ⇡ r k p(y w,i ; r k ) p(y w,i ; # r ) .
M step: Updating the mixture parameter estimates by maximizing the expected complete log-likelihood

⇡ r +1 k = P i t ik (# r ) n and ( jk ) r +1 = P n i=1 t ik (# r )y jẁ ,i P n i=1 t ik (# r )
.

The stopping rule can rely on a plateau of the log-likelihood (equally L(#; yw) or L(#|w; x)) or a given iteration number. Several runs should be performed also from di↵erent starting parameters for avoiding local maxima traps. is the ICL criterion for LCM. This latter is classically defined as follows

ICL LCM w = K X k=1 d X j=1 8 < : L j X `=1 ln ⇣ njẁ ,k + 1 2 ⌘ ln (n w ,k + L j 2 ) 9 = ; ln (n + K 2 ) + ln ( K 2 ) +K d X j=1 n ln ( L j 2 ) L j ln ( 1 2 ) o + K X k=1 ln (n w ,k + 1 2 ) K ln ( 1 2 ),
where nw,k = #{i : ẑw,ik = 1} and njẁ ,k = #{i : ẑw,ik = 1, y jẁ ,i = 1}.

As a consequence, the ICL LCM-LM w is itself easily calculated also.

  Classical model-based clustering for HD-d methodsThe need of model-based clustering for HD-m j methods5/21The general trade-o↵ bias/variance as a guidelineThe fundamental statistical principle Always minimize an error err between truth (z) and estimate (ẑ) Gap between true (z) and model-based (Zp) partitions: z ⇤ = arg min z2Z p err(z, z) Estimation ẑ of z ⇤ in Zp: any relevant method (bias, consistency, e ciency. . . ) Fundamental decomposition of the observed error err(z, ẑ): err(z, ẑ) = n err(z, z ⇤ ) err(z, z) o + n err(z, ẑ) err(z, z ⇤ ) Example for HD-d: reduce variance, accept bias A two-component d-variate Gaussian mixture with intra-dependency:

  Other classical ways for reducing variance in HD-d Dimension reduction In non-canonical space (PCA-like typically): Gaussian mixture of factor analysers [Ghahramani and Hinton, 97], [McLachlan et al., 03], [McNicholas and Murphy, 10], HD Gaussian models [Bouveyron et al., 07] In the canonical space (variable selection): Gaussian "variable selection" [Raftery and Dean, 06], [Maugis et al., 09a], [Maugis et al., 09b], [Sedki et al., 14]Model parsimony (constraints on model parameters): previous example, co-clustering[Govaert, 2011],[Biernacki et al., 23]

  : a LCM extension for the case of numerous levels Idea: constrain some levels to share the same parameter value (level partitionning) Notation: for each variable j, L j clusters of levels w j = (w 1 j , . . . , w m j j ) New model LCM-LM (Latent Class Model by Levels Merging): from the full level clustering w = (w 1 , . . . , w d ), the new LCM is defined by p(x i |w; #) = 1 , . . . , K ), k = ( jk ; j = 1, . . . , d; `= 1, . . . , L j ) and # = (⇡, )An interesting extreme case: if L j = 1 then LCM-LM acts as a variable selection (no e↵ect of variable j on the clustering) LCM-LM is parsimonious with regards to the number of levels since its number of parameters is equal to K 1 the log-likelihood of LCM-LM as a classical LCM. . . LCM-LM , LCM with the following new data yw = (y w,1 , . . . , yw,n), and where y w,i = (y jẁ ,i ; j = 1, . . . , d; `= 1, . . . , L j ) with y jẁ ,i = Pm j h=1 w hj x jh i Parameter estimation: L(#|w; x) = L(#; yw) + cw , where cw = with regards to the parameter # Since cw is independent of # then use a classical LCM EM

  ICL criterion[Biernacki et al., 00],[Biernacki et al., 11]: just written as a penalization of the classical LCM ICL ICL LCM-LM w = ln p(x, ẑw|w) = ICL LCM w + cw |{z} ! where ICL LCM w

  Theoretical error decreases when d grows: err theo = ( kµ 2 µ 1 k ⌃ 1 /2)

	Empirical error rate with the (true) intra-correlated model worse with d
	Empirical error rate with the (false) intra-independent model better with d!
		4											0.38
		3											0.36
		2											0.34
																					Empirical corr.
	x2	1										err	0.32								Empirical indep. Theoretical
		0											0.3
		-1											0.28
		-2											0.26
		-4 -3	-3	-2	-1	0	x1	1	2	3	4	5	0.24	1	2	3	4	5	d	6	7	8	9	10

  LCM-LM proposal: since level 4 of variable PF is under-represented (only 2 individuals), we propose to merge it with level 3 of the same variable Results with K = 2: great improvement both in ICL and error values

	HD data clustering	Our proposal: LCM-LM Our proposal: LCM-LM	Early numerical experiments Early numerical experiments	Concluding remarks Concluding remarks
		Outline Prostate cancer data: LCM-LM vs LCM Outline Conclusion
	EM-like package: RMixtcomp 2 for LCM categorical data and also dealing with
	missing values
	Levels' distribution: we see below the frequency of each level for all these
	HD data clustering: explored and less explored cases variables (missing values are denoted by a "?" and appear below as a level, but mathematical speaking it is not a level obviously) HD data clustering: explored and less explored cases
	PF	HX	EKG	BM
	?: 4 Very promizing preliminary numerical results ?: 4 1 :168 ?: 4
	Our proposal: Latent Class Model by Levels Merging (LCM-LM) 1:450 1:289 5 :150 1:420 2: 37 2:213 6 : 75 2: 82 Our proposal: Latent Class Model by Levels Merging (LCM-LM) Need to still implement MICL
	Early numerical experiments 3: 13 3 : 51 4: 2 4 Extend simultaneously to high dimension both for the number of variables and : 26 2 (Other): 13 Early numerical experiments :2 3 the number of levels (merge co-clustering and LCM-LM?)
	Concluding remarks Concluding remarks
				LCM	LCM-LM
			ICL	-1728.4	-1618.4
			error	49.2%	28.8%
		20/21			21/21	15/21

https://cran.r-project.org/web/packages/RMixtComp/index.html 14/21

https://cran.r-project.org/web/packages/RMixtComp/index.html 19/21

Model selection: MICL criterion

MICL criterion (Maximum Integrated Complete-data Likelihood) [Marbac and Sedki,15] is a variant of ICL avoiding multiple parameter estimation

Then, the model w ⇤ maximizing MICL LCM-LM w is retained: w ⇤ = arg max w2W MICLw. Finally a specific two step algorithm derived from [Marbac and Sedki,15] can be used for performing this optimization (not implemented yet).

We do not detail this algorithm here but the idea is the following. Starting from a value w 0 uniformly sampled in the corresponding space and then a value z 0 being deduced from the MAP rule of the associated maximum likelihood estimate, iteration s > 0 of the algorithm is composed by the following two steps: Partition step Fix z s+1 such that ln p(x, z s+1 |w s ) ln p(x, z s |w s ).

Model step Fix w s+1 such that ln p(x, z s+1 |w s+1 ) ln p(x, z s+1 |w s ).

16/21

Prostate cancer data: description Individuals: 506 patients on the basis of petrial variates alone for the prostate cancer clinical trial data of [Byar and Green,80], described by four categorical variables with various numbers of levels: performance rating (so-called PF, with 4 levels), cardiovascular disease history (HX, 2 levels), electrocardiogram code (EKG, 7 levels) and bone metastases (BM, 2 levels). There exists two stages of the desease (Stage 3 and Stage 4). Some missing data: 62 missing values (⇡ 1%)

We forget the classes (Stages of the desease) for performing clustering