# A new approach in mixed distributions detection

Edwin Diday, Anne Schroëder

► To cite this version:

Edwin Diday, Anne Schroëder. A new approach in mixed distributions detection. [Research Report] IRIA-RR-052, IRIA. 1974, pp.38. hal-04715934

HAL Id: hal-04715934

https://inria.hal.science/hal-04715934v1

Submitted on 1 Oct 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
RIA
laboria

laboratoire de recherche
en informatique
et automatique

A new approach
in mixed distributions
detection

E.Diday, A.Schroeder

38 P

Institut de Recherche
d'Informatique
et d'Automatique

Domaine de Voluceau
Rocquencourt
78150 - Le Chesnay
Tél. 954.90.20

Rapport de Recherche n° 52

Janvier 1974
# A NEW APPROACH IN MIXED DISTRIBUTIONS DETECTION

E. Diday, A. Schroeder

## Résumé

On se propose de détecter s'il existe, dans une population multidimensionnelle donnée, des sous-populations issues de distributions de type connu.

Dans une première partie, la méthode des Nuées Dynamiques est décrite dans le cas général de noyaux quelconques.

Dans une deuxième partie, elle est appliquée au problème cité plus haut, de reconnaissance de mélanges, en prenant des lois de probabilité comme noyaux.

Tel qu'il est présenté ici, l'algorithme maximise un critère de vraisemblance.

## Abstract

Our aim is to detect in a given multivariate sample the possible presence of sub-samples drawn from probability distributions of some known type.

In the first paragraph, the Dynamic Clusters method is described in the general case of unspecified kernels. Afterwards, it is applied to the above quoted problem of mixed distributions detection, taking probability distribution functions as kernels.

As presented here, the algorithm maximizes a likelihood criterion function.
TABLE OF CONTENTS

# 1 - INTRODUCTION

1.1. The problem

1.2. Various approaches

1.3. The Dynamic Clusters approach

# 2 - THE DYNAMIC CLUSTERS (D.C.) ALGORITHM

2.1. Notation

2.2. The algorithm

# 3 - MIXED DISTRIBUTIONS DETECTION

3.1. The problem

3.2. The algorithm

3.3. Convergence of the algorithm

# 4 - PRACTICAL ASPECTS AND INTERPRETATION

4.1. Meaning of the k classes obtained

4.2. The overall distribution and the resolution of mixtures

4.3. The initial partition

4.4. Notion of strong pattern. The number k.

# 5 - THE PARTICULAR CASE OF GAUSSIAN DISTRIBUTIONS

5.1. Definition of D,R,W, f and g in the Gaussian case

5.2. Geometrical interpretation

5.3. Sebestyen's problem

# 6 - FOUR NUMERICAL EXPERIMENTS

# 7 - CONCLUSION

## REFERENCES
# 1 - INTRODUCTION

## 1.1. - THE PROBLEM

From a multivariate sample, all multidimensional techniques (principal components, factor analysis, multidimensional scaling, clustering techniques, hierarchical or not) give a description of the population, that have afterwards to be interpreted either by further statistical methods (e.g. discrimination,...) or by common sense. These techniques need no probabilistic assumption, neither necessarily refer to any probabilistic model.

On the other hand, problems and techniques of modelization exist which try to adapt a stochastic model to a real phenomenon.

Our purpose can be situated between the two above : it consists in detecting if a given sample contains sub-samples which could be drawn from populations distributed according to a known family of probability distribution functions.

More precisely, given an empirical distribution (for instance, a histogram in the univariate case), we shall try to find if it results from the "effects" of several stochastic phenomena arising from different distributions.

For instance, the problem is the classical Resolution of Mixtures when one can assume that the overall distribution function - which will be written F(x) - actually follows the model :

$$ (1) \quad F(x) = \sum_{j<k} P_j \cdot F_j(x) $$

- $P_j$ : a priori probability of the j-th distribution

- $F_j(x)$ : j-th distribution function belonging to a known family
In terms of signal theory, this same problem can be expressed (in the case : k = 2) as "a model of observation signal consisting of a mixture of two unknown pulse waveforms of some duration T, which occur independently successively at random with probabilities P₁ and P₂, imbedded in additive zero-mean stationary Gaussian noise with unknown power" - (from W.D. GREGG and J.C. HANCOCK (1969)).

## 1.2. - VARIOUS APPROACHES

A number of techniques have been proposed to this end. We shall classify them in two categories, according to the problem they solve.

(a) Those, that use model (1) knowing the number k and the form of the $F_j$, and that estimate the $P_j$ and the unknown parameters on which the $F_j$ depend :

$$F_j(x) = \varphi(x, \theta_j), \quad \theta_j \in \mathbb{R}^s \\ \theta_j \text{ unknown, but } \varphi \text{ known}$$

(b) Those that search the same model (1) and estimate k, the expressions of the $F_j$ and the $P_j$.

The first category techniques are estimation techniques, the second category includes bayesian methods, supervised or unsupervised learning, stochastic approximation, etc...

The first category techniques are estimation techniques and they only differ in the type of estimators they use : method of moments (Pearson (1894)) with maximum likelihood estimates (RAO (1948), DAY (1969)), minimum $\chi^2$, etc... Most of them require Gaussian distributions and often the univariate case ; let us note two of them : RAO (1948), who needs two classes ; BATTACHARYA (1967), who gives a graphical method to determine the number of classes, but needs a great number of observations to be collected and the different distributions to be adequately separated ; for a general review of these techniques, see DOROFEYUK (1971).
DAY (1969) also deals with two components, that may be multi-variate.

Another approach to the estimation problem of model (1) is COOPER and COOPER'S (1964) : the unknown parameters are deduced from the moments of the overall observed distribution.

Many hypotheses must be made when studying general cases such as: when the number of components is greater than two, or for multivariate Gaussian distributions, etc... (e.g. DAY assumes that the covariance matrices are equal, from a practical point of view, COOPER and COOPER study the case of two distributions only differing by their means).

The algorithms of approximation type are very different. Most of them attempt a bayesian approach (PATRICK and HANCOCK (1966), PATRICK and COSTELLO (1970), PATRICK (1972), AGRAWALA (1970)) though giving quite various techniques. The hypotheses to be made differ from one method to another, but are usually very restrictive.

This sort of approach permits one to formalize the mixed distributions detection problem in terms of unsupervised learning (AGRAWALA (1970) PATRICK (1972), DUDA and HART (1973)).

An information criterion can also be used for a stochastic approximation algorithm (YOUNG and CORALUPPI (1970)) which is very interesting, because it does not need to know the actual number of components in the mixture ; it is however restricted to one-dimensional Gaussian distributions.

### 1.3. - THE DYNAMIC CLUSTERS APPROACH

Our own approach does not refer in any way to model (1). We only want to detect in the population the possible presence of samples of some known distribution, but we do not make any assumptions about the global distribution.

However, if the user can make hypotheses on the representative-ness of the population as a sample of a global population and can admit model (1) for the overall distribution, then our algorithm will give him a solution of the mixture type.
We shall use an algorithm of the "Dynamic clusters" type (cf. DIDAY, [6] to [9]) i.e. an algorithm that detects parallely clusters among the observations and some typical features for these clusters : the "kernels"; for us the feature to find will be fitness with a probability distribution of some known type.

Our algorithm will be presented in a general form for any distribution function, any dimension for the sample and any number of components.

The only input we need is the form of the probability distribution function and the number of components. However, we can see that this last hypothesis is not really restrictive, because the actual number can be found by several means even if the algorithm has run with another a priori number.

D.C. algorithm has been thoroughly studied in E. DIDAY ([6] to [9]) in two particular cases : when the kernels are linear manifolds (in factor analysis), and when they are subsets of the population, in non-hierarchical clustering ; in that case, the algorithm belongs to the class of methods such as Isodata (Hall and Ball - 1965), k-means (Mac Queen - 1967), iterative relocation (Wishart - 1971).

In the first part of the paper, we shall recall the D.C. (Dynamic Clusters) algorithm in the general case - for any set of kernels - and its convergence properties. Afterwards, the problem will be formalized in precise terms and DC applied to it.

Several practical remarks will then be made, essentially concerning the number of components. Afterwards, the particular case of mixed gaussian distributions will be developed, since it brings interesting results in a geometrical viewpoint. We shall then conclude with several examples.
## 2 - THE DYNAMIC CLUSTERS ALGORITHM

### 2.1 - NOTATION

- Let $E$ be any finite subset of $\mathbb{R}^q$.

- $P_k$: the set of all partitions of $E$ into $k$ classes;

the elements of $P_k$ will be called : k-partitions.

$$PeP_k \iff P = (P_1, \dots, P_k)$$

- $\mathcal{L}$: a set that will be called the space of "kernels". These kernels will be associated with subsets of $E$, as a characterization of these subsets depending on the application of the algorithm.

- $\mathcal{L}_k$: the set of all k-tuples of $\mathcal{L}$:

$$ \begin{gathered} L \in \mathcal{L}_k \iff L = (\lambda_1, \lambda_2, \dots, \lambda_k) \quad \text{where } \lambda_i \in \mathcal{L} \quad \forall i \in k] \\ (\text{if } k] = \{1, 2, \dots, k\}) \end{gathered} $$

The D.C. algorithm solves the following problem : To find a couple $(P,L)$, where $PeP_k$ and $L \in \mathcal{L}_k$ that minimizes some criterion function, which will be denoted $W : \mathcal{L}_k \times P_k \to \mathbb{R}^+$. Many problems can be written in these terms as long as any choice is possible for the kernels and the criterion function.

The general idea of the algorithm is then quite simple :
It consists in deducing from any element of $\mathcal{L}_k$ an element of $P_k$,
and from this element of $P_k$ an element in $\mathcal{L}_k$, so that the values of
the criterion function on successive couples $(L,P)$ decrease. To be
precise, we shall need some more notation :

$$ \bullet \ D : E \times \mathcal{L} \to \mathbb{R}^+ $$

$(x, \lambda) \leftrightarrow D(x, \lambda)$ will express a "distance" between an observation and a kernel.
$$ \bullet \ R : \mathbb{L} \times ]k] \times \mathbb{P}_k \to \mathbb{R}^+ $$

($\lambda, i, P$) $\to$ $R(\lambda, i, P)$, will give a measure of the fitness of kernel $\lambda$ with the ith element $P_i$ of the k-partition P.

The way to associate a k-partition with a k-tuple of kernels will then be given by the following function :

$$ \bullet \ f : \mathbb{L}_k \to \mathbb{P}_k $$

$$ L = (\lambda_1, \dots, \lambda_k) \implies f(L) = P = (P_1, \dots, P_k) $$

where :

$$ P_i = \{x \in E / D(x, \lambda_i) < D(x, \lambda_j), \forall j \in ]k]\} $$

(In case of equality, x will be assigned to the lower index class).

$P_i$ is therefore built with all elements of E that are "nearer" (in D sense) to $\lambda_i$ than to any other kernel of L.

Reciprocaly, to associate a k-tuple of kernels to a k-partition, we shall introduce the function g :

$$ \bullet \ g : \mathbb{P}_k \to \mathbb{L}_k $$

$$ P = (P_1 \dots P_k) \implies g(F) = L = (\lambda_1 \dots \lambda_k) $$

where, $\forall i \in ]k]$, $\lambda_i$ is given by :

$$ R(\lambda_i, i, P) = \min_{\lambda \in \mathbb{L}} R(\lambda, i, P) $$

In other words, among all possible kernels, $\lambda_i$ is chosen as the nearest (in terms of $R$) to $P_i$. (If this definition leads to several L, one has to define a unique choice, so that, for any $P \in \mathbb{P}_k$, $g(P)$ is a well determined element of $\mathbb{L}_k$).

The criterion to minimize will be defined as :

$$ \bullet \ W : \mathbb{L}_k \times \mathbb{P}_k \to \mathbb{R}^+ $$

$$ v = (L, P) \implies W(v) = \sum_{i \in ]k]} R(\lambda_i, i, P) $$
The problem is now an optimization problem :

(2) Minimize $W(v)$ for all $v \in \mathcal{L}_k \times \mathcal{P}_k$

## 2.2 - THE ALGORITHM

The D.C. algorithm is based on two sequences :

$$- (v_n), \text{ in } \mathcal{L}_k \times \mathcal{P}_k, \text{ i.e. } v_n = (L^{(n)}, P^{(n)})$$

$$- \text{ and } u_n = W(v_n) \in \mathbb{R}^+$$

Let $F^{(\circ)}$, be any initial $k$-partition (it can be drawn at random, either chosen) and $L^{(\circ)} = g(P^{(\circ)}) \rightarrow v_0 = (L^{(\circ)}, P^{(\circ)})$ the sequence $(v_n)$ is then defined recursively :

$v_{n+1} = (L^{(n+1)}, P^{(n+1)})$ is deduced from $v_n$ by : $P^{(n+1)} = f(L^{(n)})$

and $L^{(n+1)} = g(P^{(n+1)})$.

We shall show that - under certain constraints - the sequence $u_n = W(v_n)$ decreases. As it is a sequence in $\mathbb{R}^+$, it converges, and we shall see (theorem 2) that its limit is attained :

$$\exists M : Vn \ge M \quad u_n = u^*$$

A couple $v^* = (L^*, P^*)$ such as $W(v^*) = u^*$ will be called a LOCAL OPTIMUM for the problem.

A couple $(L^*, P^*)$ will be called a GLOBAL OPTIMUM if :

$$W(v^*) \le W(v) \text{ for all } v \in \mathcal{L}_k \times \mathcal{P}_k,$$

For a $v^*$ given by the algorithm, this inequality only holds for $v$ in a part of $\mathcal{L}_k \times \mathcal{P}_k$, this is why it is called a "local" optimum. For further information on this optimality, see DIDAY [9].
We shall see now how $(u_n)$ decreases.

**Definition 1:** A function $R: [L \times ]k] \times P_k \to R^+$ is said **SEMI-SQUARE**, if :

(1) $\forall L \in L_k \quad W(g \circ f(L), f \circ g \circ f(L)) \le W(g \circ f(L), f(L))$ (The sign $\le$ is used to denote the composition of functions).

(1) can be written equivalently :

$$ (1') \quad \forall P \in L_k \quad W(g(P), f \circ g(P)) \le W(g(P), P) $$

**Theorem 1:** If $R$ is semi-square, the sequence $(u_n)$ decreases and therefore converges.

**Proof:** We shall prove the following inequality :

$$ u_{n+1} \le u_n \quad \forall n $$

in two steps : a/ $u_{n+1} = W(L^{(n+1)}, P^{(n+1)}) \le W(L^{(n)}, P^{(n+1)})$

and b/ $W(L^{(n)}, P^{(n+1)}) \le W(L^{(n)}, P^{(n)}) = u_n$

a/ $u_{n+1} = \sum_{i \in ]k]} R(\lambda_i^{(n+1)}, i, P^{(n+1)})$

and $W(L^{(n)}, P^{(n+1)}) = \sum_{i \in ]k]} R(\lambda_i^{(n)}, i, P^{(n+1)})$

Since, for each $i \in ]k]$, $\lambda_i^{(n+1)}$ is deduced from $P^{(n+1)}$ by application $g$, so that :

$$ R(\lambda_i^{(n+1)}, i, P^{(n+1)}) = \min_{\lambda \in L} R(\lambda, i, P^{(n+1)}) $$

we have in particular :

$$ \forall i \in ]k] \quad R(\lambda_i^{(n+1)}, i, P^{(n+1)}) \le R(\lambda_i^{(n)}, i, P^{(n+1)}) $$

and therefore, a/ is proved.

b/ As soon as $n \ge 1$, $P^{(n)} = f(L^{(n-1)})$ and then $P^{(n)} \in f(L_k)$. Since $R$ is assumed to be semi-square, we can apply the property (1') of definition 1 to $P^{(n)}$ :
$$W(g(P^{(n)}), f \circ g(P^{(n)})) \leq W(g(P^{(n)}), P^{(n)})$$

$$\text{and as : } g(P^{(n)}) = L^{(n)} \text{ and } f \circ g(P^{(n)}) = f(L^{(n)}) = P^{(n+1)}$$

$$\text{we have : } W(L^{(n)}, P^{(n+1)}) \leq W(L^{(n)}, P^{(n)}), \text{ which proves b/}$$

Then, $(u_n)$ decreases and as it is a sequence in $\mathbb{R}^+$, it converges.

(qed)

**Remark :** The property of $\mathbb{R}$ to be semi-square is necessary and sufficient to prove b/, but it is only sufficient to have $(u_n)$ decreasing since we could have $u_{n+1} \leq u_n$ without the intermediate inequalities we have used.

**Proposition :** The two following conditions (2) and (3) on $\mathbb{R}$ are sufficient to have $\mathbb{R}$ semi-square. Moreover (2) implies (3).

$$\text{(2)} \quad \left. \begin{array}{l} V L \in L_k \\ V P \in P_k \end{array} \right\} W(L, f(L)) \leq W(L, P)$$

$$\text{(3)} \quad V L, M \in L_k : W(L, f(M)) \leq W(M, f(M)) \implies W(L, f(L)) \leq W(L, f(M))$$

If (3) is true $\mathbb{R}$ is said to be **SQUARE**.

**Proof :** One can see easily that (2) $\implies$ (1) and (2) $\implies$ (3)

Let us prove that (3) $\implies$ (1) :

As (3) stands for all $L$ and $M$, we can write it for $L = g \circ f(M)$, then :

$$\text{(3)} \quad \implies \quad V M \in L_k :$$

$$W(g \circ f(M), f \circ g \circ f(M)) \leq W(M, f(M)) \implies W(g \circ f(M), f \circ g \circ f(M)) \leq W(g \circ f(M), f(M))$$

By definition of $g$ and $W$ from $\mathbb{R}$, we have :

$$V M \in L_k, V P \in P_k \quad W(g(P), P) \leq W(M, P)$$
(This have in fact been proved in the inequality a/ of theorem 1).

Therefore, the first member of the implication is always true
and the implication is reduced to its second member, for all M, which
is exactly the condition (1) - qed -

**Corollary 1 :** If R is defined as : ∀λ ∈ L and ∀P ∈ P<sub>k</sub>,

R(λ, i.P) = Σ_{x∈P<sub>j</sub>} D(x,λ), then the sequence (u<sub>n</sub>) decreases and converges.

**Proof :** We shall see that the condition (2) is true for such an R :

Let us take L = (λ₁...λₖ) and f(L) = Q = (Q₁,...,Qₖ)

Then, for all P ∈ Pₖ :

$$
-W(L,P) = \sum_{j \in ]k]} \sum_{\{x \in P_j\}} D(x, \lambda_j) = \sum_{\{x \in E\}} \sum_{j \in ]k]} D(x, \lambda_j) \cdot \delta_{P_j}(x)
$$

(where $\delta_{P_j}(\cdot)$ is the characteristic function of $P_j$, i.e :

$$
\delta_{P_j}(x) = 1 \text{ if } x \in P_j \text{ and } = 0 \text{ otherwise}
$$

On the other hand :

$$
\begin{align*}
- W(L, f(L)) &= \sum_{i \in ]k]} R(\lambda_i, i, Q) = \sum_{i \in ]k]} \sum_{\{x \in Q_i\}} D(x, \lambda_i) \\
&= \sum_{\{x \in E\}} \sum_{i \in ]k]} D(x, \lambda_i) \cdot \delta_{Q_i}(x)
\end{align*}
$$

As Q = f(L), and by definition of f, we know that x∈Qᵢ if

$$
D(x, \lambda_i) \le D(x, \lambda_j)
$$

for all j, which implies, for all L ∈ L<sub>k</sub> and P ∈ P<sub>k</sub> : W(L,f(L)) ≤ W(L,P)

- qed -

As we have studied the convergence of (u_n), we shall now see
how (v_n) can give a solution to the problem (2).

**Definition** : An element v = (L,P) ∈ L_k x P_k is said to be **UNBIASED** for the functions f and g, if P = f.o.g(P) and L = g.o.f(L).
**Definition** : A sequence $(L^{(n)}, P^{(n)})$ in $L_k \times P_k$ is said to be convergent if there exists a $M$ such that :

$$ \forall n \ge M (L^{(n)}, P^{(n)}) = (L^*, P^*) $$

**Theorem 2** : If $R$ has the two following properties : (i) semi-square, (ii) $\forall P \in (L_k), \forall i \in P$, $R(\lambda, i, P)$ is minimum for a unique $\lambda_o$. Then, $(v_n)$ is convergent and its limit is an unbiased element.

**Proof** : Since $E$ is assumed to be a finite set, $P_k$ is finite too. Since by definition, $E$ is such that, for any $P \in P_k$, $g(P)$ can only take one value in $L_k$, then $g(P_k) \subseteq L_k$ is finite too and $(v_n)$ and $(u_n)$ can only take a finite number of values. As $(u_n)$ converges, its limits $u^*$ is attained :

$$ \exists M : \forall n \ge M \quad u_n = u^* $$

Thus, $\forall n \ge M$ $W(L^{(n)}, P^{(n)}) = W(L^{(n+1)}, P^{(n+1)})$ which implies that the two inequalities a/ and b/ of theorem 1 are equalities :

$$ W(L^{(n)}, P^{(n)}) = W(L^{(n)}, P^{(n+1)}) = W(L^{(n+1)}, P^{(n+1)}) $$

The second equality may be written :

$$ \sum_{i \in ]k]} R(\lambda_i^{(n)}, i, P^{(n+1)}) = \sum_{i \in ]k]} R(\lambda_i^{(n+1)}, i, P^{(n+1)}) $$

The hypothesis (ii) implies :

$$ \forall i \in ]k] \quad R(\lambda_i^{(n+1)}, i, P^{(n+1)}) \le R(\lambda_i^{(n)}, i, P^{(n+1)}) $$

with the equality if and only if $\lambda_i^{(n+1)} = \lambda_i^{(n)}$.

We then have two sums of positive numbers which are equal while every term of one of them is lower than or equal to the corresponding term of the other. This is only possible if all corresponding terms are equal and therefore :

$$ \lambda_i^{(n+1)} = \lambda_i^{(n)} \quad \forall i \in ]k] \Leftrightarrow L^{(n+1)} = L^{(n)} $$

Then $P^{(n+1)} = f(L^{(n)})$ and $P^{(n+2)} = f(L^{(n)})$ imply that $P^{(n+1)} = P^{(n+2)}$ and the convergence of $(v_n)$ is proved :

$$ \forall n \ge M+1 \quad v_n = v_{n+1} = v^* $$
$$v^* = (L^*, P^*) = (L^{(n)}, P^{(n)}) = (g(P^{(n)}), P^{(n)}) \implies \begin{cases} P^* = P^{(n)} \\ L^* = g(P^{(n)}) = g(P^*) \end{cases}$$

and

$$v^* = (L^*, P^*) = (L^{(n+1)}, P^{(n+1)}) = (g \circ f(L^{(n)}), f(L^{(n)})) \implies \begin{cases} P^* = f(L^*) = f \circ g(P^*) \\ L^* = g \circ f(L^*) \end{cases}$$

Therefore $v^*$ is an unbiased element for $f$ and $g$.

q.e.d.
## 3 - MIXED DISTRIBUTIONS DETECTION

### 3.1 - THE PROBLEM

We shall now write in mathematical terms the problem we have informally described in the introduction.

Let $E$ be a set of $N$ observations on which $q$ measures have been taken : then $E$ is a finite subset of $\mathbb{R}^q$.

Suppose we are given a family of distribution functions :

$$ (f_{\lambda})_{\lambda \in \mathcal{L}} $$

which depends on the parameter $\lambda$, with $\lambda \in \mathcal{LCR}^s$. (For instance, if $q = 1$, this family could be that of Gaussian univariate distributions with

$$ \lambda = (\mu, \sigma), \ s = 2, \ \mathcal{L} = \mathbb{R} \times \mathbb{R}^+ $$

We want then to find a couple $(L, P)$, where $L = (\lambda_1...\lambda_k)$, $\lambda_i \in \mathcal{L}$ and $P = (P_1...P_k)$ is a $k$-partition of $E$, such that for all $i \in \mathcal{J}$, $P_i$ may be considered as a "likely" sample of the distribution $f_{\lambda_i}$:

To this end, we shall try to maximize the product of the likelihoods of the $k$ "samples" $P_i$ for the distributions $f_{\lambda_i}$, or, in other words, to find $L^*$ and $P^*$ such that:

$$ \mathcal{L}(L^*, P^*) = \max_{\substack{L \in \mathcal{L}_k \\ P \in P_k}} \prod_{i \in \mathcal{J}} \prod_{k \in \mathcal{K}} f_{\lambda_i}(x) $$

Let us now note : $V_{\lambda_i}(P_i) = \prod_{k \in \mathcal{K}} f_{\lambda_i}(x)$;

We shall show that, in fact, the D.C. algorithm make the following criterion(*) on $(L, P)$ decrease :

$$ W(L, P) = K - \sum_{i \in \mathcal{J}} \log V_{\lambda_i}(P_i) $$

(where $K$ is a constant)

<footnote>(*) The same problem may of course be formalized in other terms.</footnote>
3.2 - THE ALGORITHM

Let us take :

. $E \in \mathbb{R}^q$ - The finite set to classify

. $L \in \mathbb{R}^s$ - Set of the kernels (which will be exactly the set of parameters introduced above.)

$$
. D : E \times \mathcal{L} \rightarrow \mathbb{R}^{+}
$$

$$
(x, \lambda) \rightarrow D(x, \lambda) = \log (f^* / f_{\lambda}(x))
$$

where $f^* \ge \max\{f_\lambda(x) | \lambda \in \mathcal{L}, x \in E\}$, it is sufficient to know
such an $f$ exists to use the algorithm ;

for instance, for univariate Gaussian distributions with $\lambda = (\mu, \sigma)$:

$$
\max_{x \in \mathbb{R}} f_{\lambda}(x) = (2\Pi)^{-1/2} \sigma^{-2} \quad \text{and then a possible } f^{*} \text{ is:}
$$

$$
f = (2\Pi)^{-1/2} d_o^{-1} \quad \text{if } d_o = \begin{cases} \min_{x \in E} & |x - y|^2 \\ \min_{y \in E} & |x - y|^2 \end{cases}
$$

This definition for D expresses that the **greater** $f_\lambda(x)$ is, the
nearer to kernel $\lambda$ the observation $x$ is.

(It can also be said that the likelihood of the sample $\{x\}$ is large for
$f_\lambda$).

. R is then defined from D :

$$
\mathbb{R} : \mathcal{L} \times [k] \times \mathbb{P}_k \to \mathbb{R}^+
$$

$$
(\lambda, i, P) \rightarrow \mathbb{R}(\lambda, i, P) = \sum_{\{x \in P_i\}} D(x, \lambda)
$$

$$
\begin{align*}
\mathbb{R}(\lambda, i, P) &= \sum_{\{x \in P_i\}} \text{Log}(f^*/f_\lambda(x)) = \text{Log} \left[ \left( \frac{f^*}{\prod_{\{x \in P_i\}} f_\lambda(x)} \right) \right] \quad (\text{where} \\
&= \text{Log} \left[ \left( \frac{f^*}{|P_i|} \right) / \nu_\lambda(P_i) \right]
\end{align*}
$$
Finally, the criterion is :

$$
W : \mathcal{L}_k \times \mathcal{P}_k \rightarrow \mathbb{R}^+
$$

$$
(L, P) \quad \mapsto \quad W(L, P) = \sum_{i \in ]k]} R(\lambda_i, i, P) \\
\text{(if } L = (\lambda_1, \dots, \lambda_k)\text{)}
$$

$$
W(L,P) = \sum_{i \in ]k]} \log \left[ (f^*) \frac{|P_i|}{V_{\lambda_i}(P_i)} \right]
$$

$$
W(L,P) = \log(f^*)^N - \sum_{i \in ]k]} \log V_{\lambda_i}(P_i) \quad (\text{since } \sum_{i \in ]k]} |P_i| = N)
$$

So the above formalization leads us to a criterion which ex-
presses that we shall maximize the product of the likelihoods of the k
samples Pᵢ.

The two fundamental functions of the algorithm, f and g,
become :

$$
f : \mathcal{L}_k \to \mathcal{P}_k
$$

$$
L = (\lambda_1 \dots \lambda_k) \implies P = (P_1 \dots P_k)
$$

$$
\text{where } P_i = \{x \in E \mid D(x, \lambda_i) < D(x, \lambda_j), \forall j \neq i \}
$$

$$
P_i = \{x \in E / f_{\lambda_i}(x) \ge f_{\lambda_j}(x), \forall j \neq i \}
$$

(x being assigned to the lower index class in case of quality)

The elements of E are therefore assigned to the class to which
they more likely belong.

$$
\cdot g : \mathcal{P}_k \to \mathcal{L}_k
$$

$$
P = (P_1 \dots P_k) \implies L = (\lambda_1 \dots \lambda_k) \text{ where } \lambda_i \text{ is such that :}
$$

$$
R(\lambda_i, i, P) = \min_{\lambda \in \mathcal{L}} R(\lambda, i, P)
$$

$$
\iff R(\lambda_i, i, P) = \min_{\lambda \in \mathcal{L}} \log \left[ (f^*) \frac{|P_i|}{V_{\lambda_i}(P_i)} \right]
$$

$$
\iff V_{\lambda_i}(P_i) = \max_{\lambda \in \mathcal{L}} V_{\lambda}(P_i)
$$

$\lambda_i$ is the maximum likelihood estimator of $\lambda$,
deduced from sample $P_i$.
This definition determines uniquely $g(P)$ since, in usual conditions of regularity for distribution functions, and at least when $f_λ$ is the general family of the exponential type distributions, there exists one and only one maximum likelihood estimator for $λ$. (see for instance FOURGEAUD and FUCHS).

### 3.3 - CONVERGENCE OF THE ALGORITHM

**Theorem 3 :**

Given :

- a number of classes, $k$;

- a family of distribution functions $(f_λ)_{λ∈L}$;

the sequence $\mathcal{L}(L^{(n)}, P^{(n)})$ of the products of the k likelihoods increases and converges.

The corresponding sequence $v_n = (L^{(n)}, P^{(n)})$ converges towards an unbiased element.

**Proof :**

As $R$ is defined as in the hypotheses of the corolary 1 (see 2.2), the sequence $W(L^{(n)}, P^{(n)})$ decreases, and since :

$$W(L^{(n)}, P^{(n)}) = \text{Cst} - \sum_{i \in I_k} \text{Log } V_{\lambda_i}^{(n)}(P_i^{(n)}) = \text{Cst} - \text{Log } \mathcal{L}(L^{(n)}, P^{(n)}),$$

$\mathcal{L}(L^{(n)}, P^{(n)})$ increases and converges.

The convergence of $(v_n)$ is ensured with the help of theorem 2.

q.e.d
# 4 - PRACTICAL ASPECTS AND INTERPRETATION

## 4.1 - MEANING OF THE K CLASSES OBTAINED

Given D.C. algorithm applied to mixed distributions detection, the k classes obtained are attached to the two following constraints :
- The notion of likelihood which has been taken as a quality criterion.
- The family of distributions which is initially chosen.

Before anything else, adequation tests must be made for each class i, between the sample $P_i$ and its presumed distribution function $f_{\lambda_i}$.

## 4.2 - THE OVERALL DISTRIBUTION AND THE RESOLUTION OF MIXTURES

In the particular case where E can be supposed a representative sample arising from an underlying distribution and when model(1) is assumed, then an overall distribution on $\mathbb{R}^q$ can be deduced from E and (L,P) :

$$Vz \in \mathbb{R}^q \quad F(z) = \sum_{i \in J_k} \Pr(z \in P_i) f_{\lambda_i}(z)$$

where the probabilities $Pr(z \in P_i)$ are estimated by the frequencies : $card(P_i)/N$.

The D.C. algorithm then gives a solution to the mixtures resolution problem ; this solution must of course be verified afterwards with adequation tests. (for such an application of the D.C. algorithm, see example 4.2).
## 4.3 - THE INITIAL PARTITION

It has been seen that the algorithm needs a partition of E as a starting point. The first idea which comes to mind to build this initial partition is of course a random classification. Several trials have proved that it was not possible to do so : In fact, if the k classes $P_1^{(0)}$, $P_2^{(0)}$ and $P_k^{(0)}$ are really uniformly distributed on E, the maximum-likelihood estimators of the unknown parameters will take almost equal values ; this means that $L^{(0)} = g(P^{(0)}) = (\lambda_1^{(0)}...\lambda_k^{(0)})$ will be such that: $\lambda_1^{(0)} \neq \lambda_2^{(0)} \neq ... \neq \lambda_k^{(0)}$ and the new partition $P^{(1)} = f(L^{(0)})$ will be very confused ; the algorithm has then difficulties to converge and takes anyway a large number of iterations.

Consequently, we have adopted a particular way of choosing $P^{(0)}$: We cut up the ranges of E in all q dimensions according to the given number k.

- If $q = 2, k = 3$ and given an $a \in ]0,1]$, the cutting up is shown in figure 1.

<figure><img src="image_2.png" /><figcaption>- fig. 1 -</figcaption></figure>

$$X^M \in E \quad \text{such as :} \quad X_1^M + \alpha X_2^M = \begin{cases} \max & (x_1 + \alpha x_2) \\ a & x \in E \end{cases}$$

$$X^m \in E \quad \text{such as :} \quad X_1^m + \alpha X_2^m = \begin{cases} \min & (x_1 + \alpha x_2) \\ a & x \in E \end{cases}$$

By this mean, the space $\mathbb{R}^p$ is cut up by hyperplanes, but one can imagine many other ways to get the initial partition : cutting up, for instance, by curves, or around spheres etc...
Moreover, it is not necessary to start from a partition, the user may have rather choosing initial kernels $L^{(0)}$; $P^{(1)}$ is then taken as $f(L^{(0)})$ and the algorithm goes on as above. A way to get an initial $L^{(0)}$ is for instance to draw at random k elements of E, take them as mean vectors for the k-initial distributions and calculate the other parameters of the distributions so that the range of E in $\mathbb{R}^p$ is entirely covered. (For instance, if $p=1$ and in the case of Gaussian distributions, the means are drawn at random and the standard deviations may be chosen so that any point of E belongs to a confidence interval of two standard-deviations around the mean for at least one of the distributions).

### 4.4 - NOTION OF STRONG PATTERN. THE NUMBER k

Knowing that different initial partitions lead to different unbiased elements, it seems necessary to compare the elements obtained from a given set of data with several trials of different initial partitions.

To this end, one has to examine the following table (fig. 2) ; which is called the "array of strong patterns" :

<figure><img src="image_3.png" /><figcaption>Fig. 2</figcaption></figure>

$n_j(x) = $ number of the class to which $x$ belongs in the partition obtained at the $j$th trial.
The elements of E are then compared with one another with the help of the following measure of dissimilarity :

$$ V_{x,y} \in E \qquad \Theta(x,y) = \begin{array}{l} \text{number of trials in which x and y have not} \\ \text{been classified in the same class.} \end{array} $$

By the way, x and y are the nearer the more often they belong to the same class of the obtained partitions. $\Theta(x,y) = 0$ means that they always belong to the same group.

**Definition :**

A strong pattern is a subset A of E, such that :

$$ V_{x,y} \in A \iff \Theta(x,y) = 0 $$

The set $\mathcal{F}$ of all strong patterns is nothing else than the quotient space E/R where R is the relation of equivalence :

$$ xRy \iff \Theta(x,y) = 0 $$

The properties of the strong patterns are thoroughly developed in DIDAY ([6] to [9]).

The example 4.4 shows how they are of use as a complement of information about the real number of classes existing in the population.
# 5 - THE PARTICULAR CASE OF GAUSSIAN DISTRIBUTIONS

## 5.1 - DEFINITION OF D,R,W, f AND g IN THE GAUSSIAN CASE

Let now the family of probability distributions ($f_{\lambda}/\lambda \in \mathcal{L}$) be the Gaussian family, i.e. :

$$V x \in \mathbb{R}^q, f_{\lambda}(x) = (2\pi)^{-q/2} (\det V)^{-1/2} \exp \left[ -\frac{1}{2} (x - \mu)^T V^{-1} (x - \mu) \right]$$

where $\lambda = (\mu, V)$ with :

* $\mu \in \mathbb{R}^q$: mean-vector of the distribution
* $V$: its covariance matrix $(q \times q)$

- Here, $\mathcal{L} = \mathbb{R}^q \times \mathcal{E}$ (if $\mathcal{E}$ is the space of all $(q \times q)$ symmetric positive definite matrices) and we have $\mathcal{L} \subset \mathbb{R}^s$ with $s = q(q+1)$.

Now, we can define D,R,W,f and g :

$$D(x, \lambda) = \text{Cst} + \frac{1}{2} \left[ \text{Log det } V + (x - \mu)' V^{-1} (x - \mu) \right]$$

- As $V^{-1}$ is a positive definite and symmetric matrix, it can be considered as a metric matrix and we shall denote by $d_V^{-1}$ the distance it defines :

$$V x, y \in \mathbb{R}^q : d_V^{-1} (x, y) = (x-y)' V^{-1} (x-y) \\ \Rightarrow D(x, \lambda) = \text{Cst} + \frac{1}{2} \left[ \text{Log det } V + d_V^{-1} (x, \mu) \right]$$

$$R(\lambda, i, P) = \sum_{\{x \in P_i\}} D(x, \lambda)$$

$$R(\lambda, i, P) = \text{Cst} + \frac{1}{2} \left[ \sum_{\{x \in P_i\}} (\text{Log det } V + d_V^{-1} (x, \mu)) \right]$$

* $\sum_{\{x \in P_i\}} d_V^{-1} (x, \mu)$ is the quadratic dispersion of the set $P_i$ around the point $\mu$ for the metric $V^{-1}$

$$W(L, P) = \text{Cst} + \sum_{i \in [k]} R(\lambda_i, i, P) \text{ where } \lambda_i = (\mu_i, V_i)$$

$$W(L, P) = \text{Cst} + \frac{1}{2} \sum_{i \in [k]} \sum_{\{x \in P_i\}} (\text{Log det } V_i + d_V^{-1} (x, \mu_i)) \\ = \text{Cst} + \frac{1}{2} \left[ \sum_{i \in [k]} |P_i| \text{Log det } V_i + \sum_{i \in [k]} \sum_{\{x \in P_i\}} d_V^{-1} (x, \mu_i) \right]$$
f : L ↦ P with :

$P_i = \{x \in E/D(x, \lambda_i) \le D(x, \lambda_j), \forall j \ne i, \text{ with } x \text{ assigned to the lower index class in case of equality } \}$

$\forall i \in ]k]: P_i = \{x \in E/\text{Log det } V_i + d_{V_i-1}^2(x, \mu_i) \le \text{Log det } V_j + d_{V_j-1}^2(x, \mu_j), \forall j\}$

g : P ↦ L, where, for all i ∈ ]k], μᵢ and Vᵢ are the maximum-likelihood estimates of the mean-vector and of the covariance matrix of the sample Pᵢ.

We know these estimates are given by :

$$ \mu_i = \frac{1}{|P_i|} \sum_{\{x \in P_i\}} x \quad \text{and} \quad V_i = \frac{1}{|P_i|} \sum_{\{x \in P_i\}} (x - \mu_i)(x - \mu_i)' $$

## 5.2 - GEOMETRICAL INTERPRETATION

In this particular case we see that the function f reclassifies the elements of E in the following way :

The "distance" between an $x \in E$ and the $i^{th}$ kernel $\lambda_i = (\mu_i, V_i)$ is expressed as the sum of two quantities :

$$ - d_{V_i-1}^2(x, \mu_i) = \text{distance from } x \text{ to } \mu_i \text{ for the metric } V_i $$

and $\text{Log det } V_i$, which does not depend on $x$ but only on $V_i$ and is a characteristic feature of the dispersion of the $i^{th}$ distribution.

Therefore $k$ kernels define on $\mathbb{R}^q$ $k$ local metrics.

The unbiased element obtained at the point of convergence gives a system of $k$ local metrics $\Delta_i$ around $k$ different points $\mu_i$ (the mean-vectors) such that :
- $\Delta_i$ is entirely defined by $\mu_i$ and a positive definite, symmetric matrix $V_i$, with the distance between $x \in \mathbb{R}^q$ and $\mu_i$ in terms of $\Delta_i$:

$$ \text{Log det } V_i + d_{V_i}^2 (x, \mu_i) $$

and:

- if $P = (P_1...P_k)$ is the k-partition of E which is determined by the $\lambda_i$, ($i = 1, ..., k$), i.e.:

$P_i = \{x \in E/x \text{ is nearer in terms of } \Delta_i \text{ to } \mu_i \text{ than to any } \mu_j \text{-in terms of } \Delta_j\}$, then, $\mu_i$ = mean-vector of $P_i$ and $V_i$ = covariance matrix of $P_i$.

We can then consider that our algorithm has given a solution to the following problem: To find local metrics in $\mathbb{R}^q$ that express in some way the features of E.

In fact, the countour-lines of the points that are equidistant from $\mu_i$ in terms of $\Delta_i$ are the ellipsoids of inertia of the Gaussian distribution of parameters $(\mu_i, V_i)$:our algorithm, in this case, is therefore able to detect ellipsoidal clusters.

## 5.3 - SEBESTYEN'S PROBLEM (see SEBESTYEN, and ROMEDER)

In his works on clustering and discrimination, Sebestyen has been brought to the following problem:

Knowing a finite population E of N elements, in $\mathbb{R}^q$, to find the distance d in $\mathbb{R}^q$ that minimizes the mean of the squares of the d-distances between all N points, two by two. If this mean is denoted by $D^2$:

$$ D^2 = \frac{1}{N(N-1)} \sum_{\{x \in E\}} \sum_{\{y \in E\}} d^2(x,y) $$
In fact, he searches d in the class of the euclidean metrics
that are defined by a positive definite symmetric matrix, and there-
fore looks for such a matrix Q which minimizes $D^2$:

$$
D^2 = \frac{1}{N(N-1)} \sum_{\{x \in E\}} \sum_{\{y \in E\}} (x - y)' Q (x - y)
$$

-NB : One remembers that any positive definite and symmetric matrix
Q can be written : Q = W'W, where W is triangular. Then, to assign
the metric Q to Rq is equivalent to assign the usual metric to
the transformed of Rq by the application : x $\rightarrow$ Wx :

$$
d_Q^2(x,y) = (x-y)'Q(x-y) = (x-y)'W'W(x-y) = (Wx-Wy)'(Wx-Wy)
$$

The problem is then solved by the following theorem :

The metric on Rq, defined by Q = W'W - where W is a linear transfor-
mation on Rq that keeps volumes constant, i.e. det W = 1 - which
minimizes D2 on the finite set E, is :

$$
Q = (\det V)^{1/q} \cdot V^{-1}
$$
where V is the covariance matrix of E.

Sebestyen's metric therefore is $d_{V^{-1}}$, with a multiplicative
constant which comes from the constraint

$Q = W'W. \text{ and } \det W = 1$

The deformation on E is then given by $V^{-1}$ : the countour-lines
of equidistant points around the mean vector $\mu$ of E are ellipsoids :

$$
\begin{align*}
d_Q^2(x, \mu) &= \text{Cst} \\
\iff (\det V)^{1/q} (x - \mu)' V^{-1} (x - \mu) &= \text{Cst} \\
\iff (x - \mu)' V^{-1} (x - \mu) &= \text{Cst}
\end{align*}
$$

which is the equation of an ellipsoid, the axis of which are given
by $V^{-1}$.
Suppose now given a k-partition P of E, Sebestyen finds local
metrics associated with the classes Pᵢ that minimize the mean of the
dispersions D²ᵢ of each class.

As for the D.C. algorithm, no k-partition is initially given.
Classes and local metrics are simultaneously researched and one can
remark that it leads to the same ellipsoidal countour lines as
Sebestyen's : the metrics differ only by constants, but not in di-
rection. (It is natural, for instance, that there is no constraint
on volumes in D.C. algorithm, since different clusters have to be
compared. It is for the same reason that it needs the additive cons-
tants : Log det Vᵢ which are associated with the dispersions of the
clusters).

On the other hand, the Sebestyen's criterion D_i^2 may be written:

if $\mu_i = \frac{1}{N_i} \sum_{\{x \in P_i\}} x$ = mean-vector of $P_i$ ($N_i = \text{Card}(P_i)$)

$$
\begin{align*}
(x - y)'Q(x - y) &= [(x - \mu_i) - (y - \mu_i)]'Q[(x - \mu_i) - (y - \mu_i)] \\
&= (x - \mu_i)'Q(x - \mu_i) + (y - \mu_i)'Q(y - \mu_i) - 2(x - \mu_i)'Q(y - \mu_i)
\end{align*}
$$

$$
\begin{align*}
D_i^2 &= \frac{2}{N_i (N_i - 1)} \left[ \sum_{\{x \in P_i\}} \sum_{\{y \in P_i\}} (x - \mu_i)' Q (x - \mu_i) - \sum_{\{x \in P_i\}} \sum_{\{y \in P_i\}} (x - \mu_i)' Q (y - \mu_i) \right] \\
&= \frac{2}{N_i (N_i - 1)} \left[ N \cdot \sum_{\{x \in P_i\}} d_Q^2 (x, \mu_i) - \sum_{\{x \in P_i\}} (x - \mu_i)' Q \left( \sum_{\{y \in P_i\}} (y - \mu_i) \right) \right]
\end{align*}
$$

and since $\sum_{\{y \in P_i\}} (y - \mu_i) = 0$, we have :

$$
D_i^2 = \frac{2}{N_i - 1} \sum_{\{x \in P_i\}} d_Q^2(x, \mu_i)
$$

Then, given a P = (P₁...Pₖ), Sebestyen minimizes the Dᵢ² by
choosing the optimal Q, while the D.C. algorithm finds simultaneously
the k-partition P and the local metrics that tend to minimize

$$
\sum_{i \in [k]} (N_i \log \det V_i + \sum_{\{x \in P_i\}} d_Q^2(x, \mu_i)) = \sum_{i \in [k]} N_i \log \det V_i + \frac{N_i - 1}{2} D_i^2
$$
This paragraph leads one to remark that the local transformations we have found in the aim of maximizing local likelihoods, while searching Gaussian distributions, belong to the family of those that minimize the mean of the square of the distances within the clusters.
**6.1** - In an univariate population drawn from three Gaussian populations, we drew :

50 observations from a population of parameters

μ₁ = 0 and σ₁ = 1,

μ₂ = 3 and σ₂ = 2,

μ₃ = -5 and σ₃ = 2,

<table><thead><tr><th>Results obtained with k=3 after<br>4 iterations</th><th>Empirical moments of the 3 drawn<br>samples</th></tr></thead><tbody><tr><td>μ₁ = -.04 σ₁ = .9</td><td>μ₁ = .1 σ₁ = 1.0</td></tr><tr><td>μ₂ = 3.5 σ₂ = 1.2</td><td>μ₂ = 3.4 σ₂ = 1.8</td></tr><tr><td>μ₃ = -4.9 σ₃ = 1.7</td><td>μ₃ = -4.9 σ₃ = 1.7</td></tr></tbody></table>

**6.2** - The data : an artificial sample proposed by DUDA and HART (1973): 25 observations drawn from the one-dimensional two components Gaussian mixture :

$$P_1 = 1/3 \qquad \mu_1 = -2 \qquad \sigma_1 = 1$$

$$P_2 = 2/3 \qquad \mu_2 = +2 \qquad \sigma_2 = 1$$

- Using the DCA algorithm, in the particular case of gaussian distributions

Input : k = 2

Output: P₁ = 8/25, μ₁ = -2.2, σ₁ = .8

P₂ = 17/25, μ₂ = 1.8, σ₂ = 1.2

(The convergence is achieved in 2 iterations, we obtained exactly the two drawn samples, associated with their maximum-likelihood parameters)

- Using Duda and Hart method (which necessarily requires Gaussian distributions)

Input : k == 2, σ₁ = σ₂ = 1, P₁ = 1/3, P₂ = 2/3

Output: μ₁ = -2.1, μ₂ = 1.7
6.3 - 150 points in R² have been drawn from three two-dimensional Gaussian populations (see fig. (a) and (b))

$$ \mu_1 = (0 ; 3) $$

$$ \mu_2 = (3 ; 0) $$

$$ \mu_3 = (3 ; 3) $$

$$ \text{and } \Sigma_1 = \Sigma_2 = \Sigma_3 = \begin{pmatrix} 1 & 0 \\ 0 & 1 \end{pmatrix} $$

<figure><img src="image_4.png" /><figcaption>Fig. (a) - Original data</figcaption></figure>

<figure><img src="image_5.png" /><figcaption>Fig. (b) - Ellipsoids of equiprobability (95 %) of the distributions from which the original data arise.</figcaption></figure>

At the point of convergence, three classes are found (see fig. (c)), the parameters of which are :

$$ \mu_1 = (-.3 ; 2.8) $$

$$ \Sigma_1 = \begin{pmatrix} .8 & -.3 \\ -.3 & .8 \end{pmatrix} $$

$$ \mu_2 = (2.9 ; -.2) $$

$$ \Sigma_2 = \begin{pmatrix} 1.1 & -.2 \\ -.2 & .8 \end{pmatrix} $$

$$ \mu_3 = (2.9 ; 3.1) $$

$$ \Sigma_3 = \begin{pmatrix} 1.3 & -.5 \\ -.5 & .9 \end{pmatrix} $$

<figure><img src="image_6.png" /><figcaption>Fig. (c) - Ellipsoids of equiprobability (95 %) of the</figcaption></figure>
This example shows the efficiency of the algorithm even if the classes are not clearly separated.

6.4 - 150 points in R² have been drawn from three two-dimensional Gaussian populations : (see fig (d) and (e)) :

$$ \mu_1 = (0; 0) \qquad \Sigma_1 = \begin{pmatrix} 4 & 2\sqrt{\frac{3}{2}} \\ 2\sqrt{\frac{3}{2}} & 1 \end{pmatrix} \qquad 2\sqrt{\frac{3}{2}} = 1.7 $$

(The principal axis of the equiprobable ellipsoids of this distribution make a $\frac{\pi}{3}$ angle with the 1st coordinate axis).

$$ \mu_2 = (0; 3) \qquad \Sigma_2 = \begin{pmatrix} 1/4 & 0 \\ 0 & 1/4 \end{pmatrix} $$

$$ \mu_3 = (4; 3) \qquad \Sigma_3 = \begin{pmatrix} 4 & -2\sqrt{\frac{3}{2}} \\ -2\sqrt{\frac{3}{2}} & 1 \end{pmatrix} \quad (\rightarrow a \frac{5}{3}\pi \text{ angle}) $$

<figure><img src="image_7.png" /><figcaption>Fig. (d) - Original data</figcaption></figure>

<figure><img src="image_8.png" /><figcaption>Fig. (e) - Ellipsoids of equiprobability (95%) of the distributions from which the original data arise</figcaption></figure>

Two trials have been done :

1/ Asking for 3 classes, the algorithm has been used starting from different initial partitions and the results achieved with the best criterion value are : (see fig (f))

$$ \mu_1 = (-.1; -.1) \qquad \Sigma_1 = \begin{pmatrix} 2.8 & 1.4 \\ 1.4 & 1.0 \end{pmatrix} $$

$$ \mu_2 = (-.1; 2.9) \qquad \Sigma_2 = \begin{pmatrix} .22 & -.03 \\ -.03 & .26 \end{pmatrix} $$
$$ \mu_3 = (3.9 ; 3.0) $$

$$ \Sigma_3 = \begin{pmatrix} 3.2 & -1.8 \\ -1.8 & 1.4 \end{pmatrix} $$

<figure><img src="image_9.png" /><figcaption>Fig. (f) - Clusters that have been detected by the algorithm after 8 iterations, asking for three classes.</figcaption></figure>

2/ Asking for 4 classes, the algorithm has been applied starting from four different initial partitions. The classification obtained at the trial that gives the best criterion value is shown on fig (g).

<figure><img src="image_10.png" /><figcaption>Fig. (g) - Clusters that have been detected by the algorithm after 8 iterations, asking for four classes</figcaption></figure>

(Two of the initial classes are recognized, the third one is cut up in two parts).

Comparing the results obtained from the four different initial drawings there are 23 strong patterns (see fig (h)) - sets of points that have been classified together in all the 4 trials, see 4.4 - but, by gathering these strong patterns as soon as they are classed together in three of the four trials, there remains only three patterns (see fig. (i)) which are the three that have been given.
<figure><img src="image_11.png" /><figcaption>Fig. (h) - "Strong patterns" got after four different drawings asking for four classes.</figcaption></figure>

<figure><img src="image_12.png" /><figcaption>Fig. (i) - Sets of points that have been classified together in three of the four drawings asking for four classes, using the array of "strong patterns".</figcaption></figure>

This example shows how the use of several initial partitions and of strong patterns can be of help when the actual number of components of the population is not known.

## 6.5 - Two possible applications to Computer systems measurement and modeling.

Although our method have not yet been tested on actual data, it may be expected that they would be applicable to the following problems proposed by E. GELENBE and J. LEROUDIER :

The problem of choosing the pages sizes in a system to minimize the fragmentation (see GELENBE E. - 1973 -). Theoretical results have been proved for a probability density function of segment size which may be expressed as a convex combination of Erlang (Γ) densities. Our algorithm could then be used to find such a combination that would be the nearest to an observed histogram.

An attempt could also be done to model the service time for different units of the system (CPU, peripherical units, etc...) and the time required by the commands issued from the users of a time-sharing system (see LEROUDIER - 1973 -). Theoretical results have not yet been proved for a distribution function of these times which is expressed as a linear combination of Erlang but such an analytic expression could be used in queuing network models.
7 - CONCLUSION

This paper introduces the general Dynamic Clusters algorithm as a useful tool in mixed distributions detection, and presents one way among many to apply it.

The interest of the proposed method can be seen in the described experiments. Anyway, this work is only a first step in that sort of application and further research has to be made on the following points :

- the initial choice of a partition or of kernels,

- the interpretation of the strong patterns ; have they a probabilistic significance in the case of distributions detection ?

- the problem of the number of classes-which is close to the precedent-,

- the choice of other criterion functions associated with other estimation techniques for the function g.

- the use of labelled samples if they are some, which would lead to a supervised learning approach.
(1) AGRAWALA A.K. - Learning with a probabilistic teacher
*IEEE Transactions on Information Theory*, vol. IT 16, n° 4
(May 1970)

(2) BALL G.H., HALL D.J. - ISODATA - A novel method of data analysis and pattern classification - Technical Report, SRI Project 5533, Stanford Research Institute, Menlo Park, California U.S.A. (1965)

(3) BHATTACHARYA C.G. - A simple method of resolution of a distribution into Gaussian components - Biometrics - (March 1967)

(4) COOPER D.B. and COOPER P.W. - Non supervised adaptive Signal detection and pattern recognition - Information and Control - 7 - pp. 416 - (1964).

(5) COOPER P.W. - Some topics on nonsupervised adaptive detection for multivariate normal distributions - Computer and Information Sciences - vol. II (1967)

(6) DAY N.E. - Estimating the components of a mixture of normal distributions - Biometrika, 56,3, pp.463 - (1969)

(7) DIDAY E. - Optimisation en classification Automatique et Reconnaissance des formes - RAIRO - Nov. 72 - V3 - pp. 61

(8) DIDAY E. - The Dynamic Clusters method in non-hierarchical clustering - International Journal of Computer and Information Sciences - (1973) - vol. 2, n° 1

(9) DIDAY E. - Une nouvelle méthode de classification automatique et reconnaissance des formes : La méthode des Nuées Dynamiques
Revue de Statistique Appliquée - (1970) - XIX - (2)

(10) DIDAY E. - Nouvelles méthodes et nouveaux concepts en classification automatique et reconnaissance des formes - Thèse de Doctorat d'Etat es-sciences Mathématiques (1972).

(11) DOROFEYUK - Automatic Classification Algorithms - (Review)
*Automation and Remote Control* - vol. 32, n° 12 part 1 (pp 1928)
(12) DUDA R.O. and HART R.E. - Pattern Classification and Scene Analysis - Wiley - New York (1973)

(13) FOURGEAUD C. and FUCHS A. - Statistiques - Dunod (Paris 1971)

(14) GELENBE E., BOEKHORST J.C.A., KESSELS J.L.W. - Minimizing wasted space in partitioned segmentation - Communications of the ACM - June 1973, vol. 16, n° 6

(15) GREGG W.D. and HANCOCK J.C. - An optimum decision - Directed scheme for Gaussian mixtures - *IEEE Transactions on Information Theory* - vol. IT-14, n° 3, (May 1968)

(16) HO Y.C. and AGRAWALA A.K. - On pattern classification algorithms; Introduction and Survey - Proceedings of the IEEE - vol. 56 n° 12 (Dec. 72)

(17) LEROUDIER J. - Une analyse de systèmes - Thèse de 3° cycle - Grenoble (Avril 1973)

(18) MAC QUEEN J. - Some Methods for Classification and Analysis of Multivariate observations - The 5th Berkley Symposium on Mathematics, Statistics and Probability - vol. 1 n° 1 (1967)

(19) PATRICK E.A. and FISCHER F.P. - Cluster mapping with experimental computer graphics - IEEE Trans. On Computers - vol c 18, N° 11 (nov. 1969)

(20) PATRICK E.A. and HANCOCK J.C. - Non supervised sequential classification and recognition of patterns - IEEE Transactions on Information Theory, vol. IT 12, n° 3 (July 66).

(21) PATRICK E.A. and COSTELLO J.P. - On unsupervised estimation algorithms - *IEEE Transactions on Information Theory*, vol. IT 16, n° 5 (sept. 1970)

(22) PATRICK E.A. - Fundamentals of Pattern Recognition - Prentice Hall INC. NJ - (1972)

(23) PEARSON K. - Contributions to the mathematic theory of evolution - *Philos. Trans. Soc.* n° 185 (1894)
(24) RAO C.R. - Utilization of multiple measurements in problems of Biological classification - *Journal of the Royal Statistical Society* - Series B - vol. X, n° 2 - (1948)

(25) ROMEDER J.M. - Méthodes de discrimination - Thèse de 3e cycle - Laboratoire de Statistique Mathématique - Université Paris VI (1969)

(26) SEBESTYEN G.S. - Decision-making process in Pattern Recognition ACM Monograph series - Macmillan, New York (1962)

(27) SPRAGINS J.S. - Learning without a teacher - *IEEE Transactions on Information Theory* - vol. IT 12, n° 2 (April 1966)

(28) WISHART D. - Some Problems in the Theory and Application of the Methods of Numerical Taxonomy - Ph. D. Thesis University of St Andrews (1971)

(29) YOUNG T.Y. and CORALUPPI G. - Stochastic estimation of a mixture of normal density functions using an information criterion - *IEEE Transactions on Information Theory*, vol. IT 16, n° 3