Model-based clustering of multivariate ordinal data relying on a stochastic binary search algorithm

We design a probability distribution for ordinal data by modeling the process generating data, which is assumed to rely only on order comparisons between categories. Contrariwise, most competitors often either forget the order information or add a non-existent distance information. The data generating process is assumed, from optimality arguments, to be a stochastic binary search algorithm in a sorted table. The resulting distribution is natively governed by two meaningful parameters (position and precision) and has very appealing properties: decrease around the mode, shape tuning from uniformity to a Dirac, identifiability. Moreover, it is easily estimated by an EM algorithm since the path in the stochastic binary search algorithm can be considered as missing values. Using then the classical latent class assumption, the previous univariate ordinal model is straightforwardly extended to model-based clustering for multivariate ordinal data. Parameters of this mixture model are estimated by an AECM algorithm. Both simulated and real data sets illustrate the great potential of this model by its ability to parsimoniously identify particularly relevant clusters which were unsuspected by some traditional competitors.

1. Ranking data occur when some judges are asking to sort m objects, quoted by \(1,\ldots ,m\), according to a preference order. Thus, a ranking data is a permutation of the whole data set \(\{1,\ldots ,m\}\). Ordinal data occur when some judges are asking to give a note among the set of ordered notes \(1<\cdots <m\). Thus, an ordinal data is a unique element of \(\{1,\ldots ,m\}\).

2. We will do the slight abuse of notation \(\hbox {p}(x;\mu ,\pi )=\hbox {p}\Big (\{x\};\mu ,\pi \Big )\).

3. For avoiding the label cluster dependency, we keep the best distance over all label permutations.

4. Note that a two component multinomial model is not performed in this univariate situation since it is not identifiable.

5. http://www.aeres-evaluation.com/

6. All partitions are better displayed in the first and the third axes than in the first and second axes.

7. We conjecture that all propositions hold for any m value. In particular, we have proved the simplest of them (Propositions (5.1), (5.2), (5.4)) in this general situation (these general proofs are not given in this paper).

Authors and Affiliations

1. Laboratoire Painlevé & Inria, Université Lille 1, 59655, Villeneuve d’Ascq, France

   Christophe Biernacki

2. Laboratoire ERIC, Université Lumière Lyon 2, 69676, Bron, France

   Julien Jacques

Corresponding author

Correspondence to Julien Jacques.

Appendix 1: The BOS properties: statements and proofs

From Equations (7), (8) and (10) in Sect. 2.2, it is easily seen that the BOS model has a polynomial expression related to \(\pi \) of degree less than or equal to \(m-1\):

From this statement, we have developed a formal Matlab code computing the exact expression of each coefficient \(a_j(m,\mu ,x)\) for any x, \(\mu \), and m values. For instance, for \(m=5\), \(\mu =2\), and \(x=4\), we obtain the following exact expression:

This formal calculus tool is used in the proofs of all the following properties. Its advantage is to provide very straightforward proofs. Its drawback is to provide proofs for only some selected values of m. However, we have chosen \(m\in \{1,\ldots ,8\}\) since it answers to most practical situations^{Footnote 7}. For instance in the reference book Agresti (2010) all ordinal data have less than 8 categories.

Proposition 5.1

(\(\pi =0\): Uniformity.) If \(\pi =0\), then \(\forall (\mu ,x) \in \{1,\ldots ,m\}^2\), \(\hbox {p}(x;\mu ,\pi )=m^{-1}\).

Proof

The formal calculus leads to \(a_0(m,\mu ,x)=m^{-1}\) for all \(m\in \{1,\ldots ,8\}\). Conclusion follows by setting \(\pi =0\) in (18).

Proposition 5.2

(\(\pi =1\): Dirac in \(\mu \).) If \(\pi =1\), then \(\forall (x,\mu ) \in \{1,\ldots ,m\}^2\), \(\hbox {p}(x;\mu ,\pi )={\mathbb I}(x=\mu )\).

Proof

The formal calculus leads to \(\sum _{j=0}^{m-1}a_j(m,\mu ,\mu )=1\) for all \(m\in \{1,\ldots ,8\}\). Conclusion follows by setting \(\pi =1\) in (18).

Proposition 5.3

(\(\mu \): mode if \(\pi >0\).) If \(\pi >0\), then \(\forall (\mu ,x) \in \{1,\ldots ,m\}^2\) such that \(\mu \ne x\), \(\hbox {p}(\mu ;\mu ,\pi ) > \hbox {p}(x;\mu ,\pi )\).

Proof

For \(\pi =1\), Proposition 5.2 already established that \(\mu \) is the mode. We consider now \(0<\pi <1\). From (18), \(\hbox {p}(\mu ;\mu ,\pi ) - \hbox {p}(x;\mu ,\pi )\) is a polynomial of degree less than or equal to \(m-1\). Moreover, from Proposition 5.1, \(\pi =0\) is a root of \(\hbox {p}(\mu ;\mu ,\pi ) - \hbox {p}(x;\mu ,\pi )\), thus we can factorize it by \(\pi \)

Since \(0<\pi < 1\) implies that \(0<\pi ^j< 1\) for any \(j>1\), the following lower bound is directly obtained

where \(b_j^-(m,\mu ,x)=\min (0,b_j(m,\mu ,x))\). All \(b_j(m,\mu ,x)\) are obtained from a formal polynomial Euclidian division for all \(m\in \{1,\ldots ,8\}\), allowing to show by formal calculus that \( b_0(m,\mu ,x) + \sum _{j=1}^{m-2} b_j^-(m,\mu ,x) > 0\) for all \(m\in \{1,\ldots ,8\}\). Conclusion follows.

Proposition 5.4

(\(\mu \): absolute growing of its probability with \(\pi \). ) \(\forall \mu \in \{1,\ldots ,m\}\), \(\hbox {p}(\mu ;\mu ,\pi )\) is an increasing function of \(\pi \).

Proof

The formal calculus leads to \(a_j(m,\mu ,\mu )\ge 0\) for all \(m\in \{1,\ldots ,8\}\) and for all \(j\in \{1,\ldots ,m-1\}\). Conclusion follows directly from (18).

Proposition 5.5

( \(\mu \): relative growing of its probability with \(\pi \).) \(\forall (\mu ,x) \in \{1,\ldots ,m\}^2\) such that \(\mu \ne x\), \(\hbox {p}(\mu ;\mu ,\pi )- \hbox {p}(x;\mu ,\pi )\) is an increasing function of \(\pi \).

Proof

The key is to prove the differential (in \(\pi \)) expression \(\hbox {p}'(\mu ;\mu ,\pi )- \hbox {p}'(x;\mu ,\pi )>0\). From (18), \(\hbox {p}'(\mu ;\mu ,\pi ) - \hbox {p}'(x;\mu ,\pi )\) is a polynomial of degree less than or equal to \(m-2\)

Since \(0<\pi < 1\) implies that \(0<\pi ^j< 1\) for any \(j>1\), the following lower bound is directly obtained

where \(b_j^-(m,\mu ,x)=\min (0,b_j(m,\mu ,x))\). All \(b_j(m,\mu ,x)\) are obtained from a formal derivation of a polynomial for all \(m\in \{1,\ldots ,8\}\), allowing to show by formal calculus again that \( b_0(m,\mu ,x) + \sum _{j=1}^{m-2} b_j^-(m,\mu ,x) > 0\) for all \(m\in \{1,\ldots ,8\}\). Conclusion follows.

Proposition 5.6

(Decreasing around \(\mu \) if \(0<\pi <1\).) \(\forall (x,x') \in \{1,\ldots ,m\}^2\) such that \(x'<x<\mu \) or \(\mu >x>x'\), \(\hbox {p}(x;\mu ,\pi )> \hbox {p}(x';\mu ,\pi )\).

Proof

From (18), \(\hbox {p}(x;\mu ,\pi ) - \hbox {p}(x';\mu ,\pi )\) is a polynomial of degree less than or equal to \(m-1\). Moreover, from Propositions 5.1 and 5.2 respectively, \(\pi =0\) and \(\pi =1\) are roots of \(\hbox {p}(x;\mu ,\pi ) - \hbox {p}(x';\mu ,\pi )\), thus we can factorize it by \(\pi (1-\pi )\)

Since \(0<\pi < 1\) implies that \(0<\pi ^j< 1\) for any \(j>1\), the following lower bound is directly obtained

where \(b_j^-(m,\mu ,x,x')=\min (0,b_j(m,\mu ,x,x'))\). All \(b_j\) \((m,\mu ,x,x')\) are obtained from a formal polynomial Euclidian division for all \(m\in \{1,\ldots ,8\}\), allowing to show by formal calculus again that \( b_0(m,\mu ,x,x') + \sum _{j=1}^{m-2} b_j^-(m,\mu ,x,\) \(x') > 0\) for all \(m\in \{1,\ldots ,8\}\). Conclusion follows.

Proposition 5.7

(Identifiability) If \(0<\pi \le 1\), identifiability holds.

Proof

The identifiability problem could concern \(\mu \) and/or \(\pi \).

• First, there exists no couple \((\mu ,\mu ')\in \{1,\ldots ,m\}^2\) with \(\mu \ne \mu '\) such that \(\hbox {p}(x;\mu ,\pi )=\hbox {p}(x;\mu ',\pi ')\) for any \(x\in \{1,\ldots ,m\}\) and any \((\pi ,\pi ')\in (0,1]^2\). Indeed, otherwise both \(\mu \) and \(\mu '\) should be modes of both distributions, which is impossible by unicity of the mode (Propositions 5.3 and 5.6).

• Second, there exists no couple \((\pi ,\pi ')\in (0,1]^2\) with \(\pi \ne \pi '\) such that \(\hbox {p}(x;\mu ,\pi )=\hbox {p}(x;\mu ,\pi ')\) for any \(x\in \{1,\ldots ,m\}\) and any \(\mu \in \{1,\ldots ,m\}\) since \(\hbox {p}(x;\mu ,\pi )\) is a polynomial of order at least one. This last statement comes from the fact that \(\hbox {p}(x;\mu ,\pi )\) varies from 1 / m to 0 or 1 when \(\pi \) varies (respectively when \(x\ne \mu \) and \(x=\mu \)) from Propositions 5.1 and 5.2, thus depends on the value of \(\pi \).

Appendix 2: AERES data

See Table 11.

Reprints and permissions