Page 70 - MSDN Magazine, November 2019

Figure 3 The Mixture Model Clustering Demo Program
The Data Structures
The demo code sets constants N = 8 (number of data items), d = 2 (data dimensionality), and K = 3 (number of clusters). There are six data structures. Matrix X is an array-of-arrays with size (8 x 2) and it holds the data. In a non-demo scenario you’d load normal- ized data into memory from a text file.
Matrix w, which has size (8 x 3), holds the membership weights and the clustering information. Each row of w represents a data item, and each column of w corresponds to a cluster. Mixture model clustering assumes that each cluster follows some probability distribution. The most commonly assumed distribution is the multi- variate Gaussian, so the technique is called Gaussian mixture model (GMM). The demo uses a simplified Gaussian, so I call the tech- nique naive Gaussian mixture model, but this isn’t a standard name.
Helper array Nk has 3 cells and holds the sums of the 3 columns of w. Array a has 3 cells and is called the mixture weights. The values in array a sum to 1.0 and are weights that indicate the contribution of each Gaussian distribution to the model. For example, if a = (0.30, 0.60, 0.10), then the first Gaussian distribution contributes 30 percent to the mixture. Each cell of the a array is initialized to 1 / 3 = 0.3333 so that each Gaussian initially contributes equally.
Matrices u and V have size (3 x 2) and hold the univariate means and variances. For example, if u[1][0] = 0.55, the mean for cluster 1, data component 0 (package height), is 0.55. If V [2][1] = 0.33, the variance for cluster 2, data component 1 (package width), is 0.33.
Updating Membership Weights
The membership weights stored in matrix w are updated in each EM iteration according to equation (3) in Figure 2. In words, for each N = 8 rows of w, you compute the K = 3 naive probabilities using the corresponding values of X, u, and V and then multiply each naive probability by the corresponding mixture weight stored in the a vector. After you compute the weighted probabilities for a row, each of the values in the row is normalized by dividing by the row sum.
Updating the membership weights is the expectation step (often shortened to E-step) of the EM optimization algorithm. After mem- bership weights in w are updated, the Nk column sums and mixture weights in vector a are updated using equation (4), the means in matrix u are updated using equation (5), and the variances in matrix V are updated using equation (6) in Figure 2. Updating Nk, a, u and V is the maximization step of EM optimization. The idea is to find the values in a, u and V that best match the data.
58 msdn magazine
Test Run
using System; namespace MixtureModel {
class MixtureModelProgram {
const int N = 8; const int K = 3; const int d = 2;
static void Main(string[] args) {
Console.WriteLine("Begin mixture model with demo "); double[][] x = new double[N][]; // 8 x 2 data
x[0] = new double[] { 0.2, 0.7 };
x[1] = new double[] { 0.1, 0.9 };
x[2] = new double[] { 0.2, 0.8 }; x[3] = new double[] { 0.4, 0.5 }; x[4] = new double[] { 0.5, 0.4 }; x[5] = new double[] { 0.9, 0.3 }; x[6] = new double[] { 0.8, 0.2 }; x[7] = new double[] { 0.7, 0.1 };
Console.WriteLine("Data (height, width): "); Console.Write("[0] "); VectorShow(x[0]); Console.Write("[1] "); VectorShow(x[1]); Console.WriteLine(". . . "); Console.Write("[7] "); VectorShow(x[7]);
Console.WriteLine("K=3, initing w, a, u, S, Nk"); double[][] w = MatrixCreate(N, K);
double[] a = new double[K] { 1.0/K, 1.0/K, 1.0/K }; double[][] u = MatrixCreate(K, d);
double[][] V = MatrixCreate(K, d, 0.01); double[] Nk = new double[K];
u[0][0] = 0.2; u[0][1] = 0.7; u[1][0] = 0.5; u[1][1] = 0.5; u[2][0] = 0.8; u[2][1] = 0.2;
Console.WriteLine("Performing 5 E-M iterations "); for (int iter = 0; iter < 5; ++iter) {
UpdateMembershipWts(w, x, u, V, a); // E step UpdateNk(Nk, w); // M steps UpdateMixtureWts(a, Nk);
UpdateMeans(u, w, x, Nk);
UpdateVariances(V, u, w, x, Nk); }
Console.WriteLine("Clustering done. \n"); ShowDataStructures(w, Nk, a, u, V);
Console.WriteLine("End mixture model demo");
Console.ReadLine(); } // Main
static void ShowDataStructures(double[][] w,
double[] Nk, double[] a, double[][] u, double[][] V)
{
Console.WriteLine("w:"); MatrixShow(w, true); Console.WriteLine("Nk:"); VectorShow(Nk, true); Console.WriteLine("a:"); VectorShow(a, true); Console.WriteLine("u:"); MatrixShow(u, true); Console.WriteLine("V:"); MatrixShow(V, true);
}
static void UpdateMembershipWts(double[][] w, double[][] x, double[][] u, double[][] V, double[] a)
{
for (int i = 0; i < N; ++i) {
double rowSum = 0.0;
for (int k = 0; k < K; ++k) {
double pdf = NaiveProb(x[i], u[k], V[k]); w[i][k] = a[k] * pdf;
rowSum += w[i][k];
}
for (int k = 0; k < K; ++k)
w[i][k] = w[i][k] / rowSum; }
}
static void UpdateNk(double[] Nk, double[][] w) {
for (int k = 0; k < K; ++k) { double sum = 0.0;
for (int i = 0; i < N; ++i)
sum += w[i][k]; Nk[k] = sum;
} }
static void UpdateMixtureWts(double[] a, double[] Nk) {
for (int k = 0; k < K; ++k) a[k] = Nk[k] / N;