EM

Expectation Maximization

In layman terms: Expectation maximization aims to find the underlying data distributions. This algorithm is a non supervised "learning" algorithm that assumes underlying Gaussian distributions and tries to iteratively find the best fitting one. As input one has the point distribution and the number of expected class centers (aka data generators). The algorithm can be helped by pushing an inital assumption to the algorithm - if nothing is provided the algorithm tries to find it's own init by randomly choosing centers and applying k-means to that distribution. Actually the difference between EM and k-means is not that big - EM simply does not assume the same covariance in all directions (a circle in 2D) but rather different sigmas in the main directions (a ellipse in 2D). The distance measure in EM is actually the Mahalonobis distance instead of simple L2.

Definition

Please refer to wikipedias great EM reference page. It contains a more mathematical definition and show some examples.

Implementation details

Creation of the object

   constructor Create( n : integer = 1000; ltol : double = 0.1);

The constructor defines the number of maximum allowed iterations to find the maximum or the under limit of error change where convergence is defined. ltol defines the percentage of the log likelihood difference between 2 iterations.

   procedure Init( W, M : IMatrix; V : IMatrixDynArr ); overload; // full init with weights, mean and covariances
   procedure Init( M : IMatrix ); overload;  // initialize centers only

These optional Init function allow to initialize the algorithm wieth either only a given initial class center M or even with initial covaricance matrices of the class centers V and the weighting matrix W. M is of size Width=x, Height=k where x is the dimension of one data point and k is the number of centers.

   function Estimate(X : IMatrix; k : integer; var W, M : IMatrix; var V : IMatrixDynArr) : boolean;

Conduct the EM algorithm on point matrix X (width=dimension of a data point. height=number of input data points) estimate weights matrix W, mean vectors M and covariance matrix V on A data point is presented row wise. k states the number of components. The function returns false if search excceeds the maximum number of iterations and throws an exception in case one center gets out of scope (aka does not get any datapoints).

If there is no init data matrix given the algorithm starts with a randomized initalization and applies a truncated k-means algorithm on the data to get a neat init.

The returned matrices contain:

• M: The centers (mean) of the found data estimation. width=number of centers, height=point dimension
• V: Array of covariance matrices of the data (length = k -> number of centers). width=height=point dimension
• W: Weighting matrix. Vector in the length of k. Values from 0 to 1 defining the number of points associated with the center.

   TExpectationMaxUpdate = procedure(Sender : TObject; iter : integer; W, M : IMatrix; V : IMatrixDynArr) of object;
   property OnUpdate : TExpectationMaxUpdate read fOnUpdate write fOnUpdate;

This event is fired on each iteration so one can visualize the algorithms performance.