Density estimation

Since the bottleneck method is framed in probabilistic rather than statistical terms, the underlying probability density at the sample points ${\displaystyle X={x_{i}}\,}$must be estimated. This is a well known problem with multiple solutions described by Silverman.^[11] In the present method, joint sample probabilities are found by use of a Markov transition matrix method and this has some mathematical synergy with the bottleneck method itself.

The arbitrarily increasing distance metric ${\displaystyle f\,}$ between all sample pairs and distance matrix is ${\displaystyle d_{i,j}=f{\Big (}{\Big |}x_{i}-x_{j}{\Big |}{\Big )}}$ . Then transition probabilities between sample pairs ${\displaystyle P_{i,j}=\exp(-\lambda d_{i,j})\,}$ for some ${\displaystyle \lambda >0\,}$must be computed. Treating samples as states, and a normalised version of ${\displaystyle P\,}$ as a Markov state transition probability matrix, the vector of probabilities of the 'states' after ${\displaystyle t\,}$ steps, conditioned on the initial state ${\displaystyle p(0)\,}$, is ${\displaystyle p(t)=P^{t}p(0)\,}$. The equilibrium probability vector ${\displaystyle p(\infty )\,}$ given, in the usual way, by the dominant eigenvector of matrix ${\displaystyle P\,}$ which is independent of the initialising vector ${\displaystyle p(0)\,}$. This Markov transition method establishes a probability at the sample points which is claimed to be proportional to the probabilities' densities there.

Other interpretations of the use of the eigenvalues of distance matrix ${\displaystyle d\,}$ are discussed in Silverman's Density Estimation for Statistics and Data Analysis.^[11]

Clusters

In the following soft clustering example, the reference vector ${\displaystyle Y\,}$ contains sample categories and the joint probability ${\displaystyle p(X,Y)\,}$ is assumed known. A soft cluster ${\displaystyle c_{k}\,}$ is defined by its probability distribution over the data samples ${\displaystyle x_{i}:\,\,\,p(c_{k}|x_{i})}$. Tishby et al. presented^[1] the following iterative set of equations to determine the clusters which are ultimately a generalization of the Blahut-Arimoto algorithm, developed in rate distortion theory. The application of this type of algorithm in neural networks appears to originate in entropy arguments arising in the application of Gibbs Distributions in deterministic annealing.^[12][13]

${\displaystyle {\begin{cases}p(c|x)=Kp(c)\exp {\Big (}-\beta \,D^{KL}{\Big [}p(y|x)\,||\,p(y|c){\Big ]}{\Big )}\\p(y|c)=\textstyle \sum _{x}p(y|x)p(c|x)p(x){\big /}p(c)\\p(c)=\textstyle \sum _{x}p(c|x)p(x)\\\end{cases}}}$

The function of each line of the iteration expands as

Line 1: This is a matrix valued set of conditional probabilities

${\displaystyle A_{i,j}=p(c_{i}|x_{j})=Kp(c_{i})\exp {\Big (}-\beta \,D^{KL}{\Big [}p(y|x_{j})\,||\,p(y|c_{i}){\Big ]}{\Big )}}$

The Kullback–Leibler divergence ${\displaystyle D^{KL}\,}$ between the ${\displaystyle Y\,}$ vectors generated by the sample data ${\displaystyle x\,}$ and those generated by its reduced information proxy ${\displaystyle c\,}$ is applied to assess the fidelity of the compressed vector with respect to the reference (or categorical) data ${\displaystyle Y\,}$ in accordance with the fundamental bottleneck equation. ${\displaystyle D^{KL}(a||b)\,}$ is the Kullback–Leibler divergence between distributions ${\displaystyle a,b\,}$

${\displaystyle D^{KL}(a||b)=\sum _{i}p(a_{i})\log {\Big (}{\frac {p(a_{i})}{p(b_{i})}}{\Big )}}$

and ${\displaystyle K\,}$ is a scalar normalization. The weighting by the negative exponent of the distance means that prior cluster probabilities are downweighted in line 1 when the Kullback–Leibler divergence is large, thus successful clusters grow in probability while unsuccessful ones decay.

Line 2: Second matrix-valued set of conditional probabilities. By definition

${\displaystyle {\begin{aligned}p(y_{i}|c_{k})&=\sum _{j}p(y_{i}|x_{j})p(x_{j}|c_{k})\\&=\sum _{j}p(y_{i}|x_{j})p(x_{j},c_{k}){\big /}p(c_{k})\\&=\sum _{j}p(y_{i}|x_{j})p(c_{k}|x_{j})p(x_{j}){\big /}p(c_{k})\\\end{aligned}}}$

where the Bayes identities ${\displaystyle p(a,b)=p(a|b)p(b)=p(b|a)p(a)\,}$ are used.

Line 3: this line finds the marginal distribution of the clusters ${\displaystyle c\,}$

${\displaystyle {\begin{aligned}p(c_{i})&=\sum _{j}p(c_{i},x_{j})&=\sum _{j}p(c_{i}|x_{j})p(x_{j})\end{aligned}}}$

This is a standard result.

Further inputs to the algorithm are the marginal sample distribution ${\displaystyle p(x)\,}$ which has already been determined by the dominant eigenvector of ${\displaystyle P\,}$ and the matrix valued Kullback–Leibler divergence function

${\displaystyle D_{i,j}^{KL}=D^{KL}{\Big [}p(y|x_{j})\,||\,p(y|c_{i}){\Big ]}{\Big )}}$

derived from the sample spacings and transition probabilities.

The matrix ${\displaystyle p(y_{i}|c_{j})\,}$ can be initialized randomly or with a reasonable guess, while matrix ${\displaystyle p(c_{i}|x_{j})\,}$ needs no prior values. Although the algorithm converges, multiple minima may exist that would need to be resolved.^[14]

An example

The following case examines clustering in a four quadrant multiplier with random inputs ${\displaystyle u,v\,}$ and two categories of output, ${\displaystyle \pm 1\,}$, generated by ${\displaystyle y=\operatorname {sign} (uv)\,}$. This function has two spatially separated clusters for each category and so demonstrates that the method can handle such distributions.

20 samples are taken, uniformly distributed on the square ${\displaystyle [-1,1]^{2}\,}$ . The number of clusters used beyond the number of categories, two in this case, has little effect on performance and the results are shown for two clusters using parameters ${\displaystyle \lambda =3,\,\beta =2.5}$.

The distance function is ${\displaystyle d_{i,j}={\Big |}x_{i}-x_{j}{\Big |}^{2}}$ where ${\displaystyle x_{i}=(u_{i},v_{i})^{T}\,}$ while the conditional distribution ${\displaystyle p(y|x)\,}$ is a 2 × 20 matrix

${\displaystyle {\begin{aligned}&Pr(y_{i}=1)=1{\text{ if }}\operatorname {sign} (u_{i}v_{i})=1\,\\&Pr(y_{i}=-1)=1{\text{ if }}\operatorname {sign} (u_{i}v_{i})=-1\,\end{aligned}}}$

and zero elsewhere.

The summation in line 2 incorporates only two values representing the training values of +1 or −1, but nevertheless works well. The figure shows the locations of the twenty samples with '0' representing Y = 1 and 'x' representing Y = −1. The contour at the unity likelihood ratio level is shown,

${\displaystyle L={\frac {\Pr(1)}{\Pr(-1)}}=1}$

as a new sample ${\displaystyle x'\,}$is scanned over the square. Theoretically the contour should align with the ${\displaystyle u=0\,}$ and ${\displaystyle v=0\,}$ coordinates but for such small sample numbers they have instead followed the spurious clusterings of the sample points.

Decision contours

Neural network/fuzzy logic analogies

This algorithm is somewhat analogous to a neural network with a single hidden layer. The internal nodes are represented by the clusters ${\displaystyle c_{j}\,}$ and the first and second layers of network weights are the conditional probabilities ${\displaystyle p(c_{j}|x_{i})\,}$ and ${\displaystyle p(y_{k}|c_{j})\,}$ respectively. However, unlike a standard neural network, the algorithm relies entirely on probabilities as inputs rather than the sample values themselves, while internal and output values are all conditional probability density distributions. Nonlinear functions are encapsulated in distance metric ${\displaystyle f(.)\,}$ (or influence functions/radial basis functions) and transition probabilities instead of sigmoid functions.

The Blahut-Arimoto three-line algorithm converges rapidly, often in tens of iterations, and by varying ${\displaystyle \beta \,}$, ${\displaystyle \lambda \,}$ and ${\displaystyle f\,}$ and the cardinality of the clusters, various levels of focus on features can be achieved.

The statistical soft clustering definition ${\displaystyle p(c_{i}|x_{j})\,}$ has some overlap with the verbal fuzzy membership concept of fuzzy logic.