Gaussian mixture modeling (Expectation-Maximization) definitions

EM initialization function

Scalar covariance initialization

EMInitializationScalar[x_, k_] := Module[{n, rl, μ, σ, priors, sm, del2}, <br />&nbs ... riors = Table[1/k, {k}] // N ; <br />     {μ, σ, priors} ] ;

Full covariance initialization

EMInitializationFull[x_, k_] := Module[{μ, σ, priors, dim},  {μ ... th[First[μ]] ;  {μ, (# IdentityMatrix[dim]) & /@ σ, priors} ] ;

Mathematica-native EM mixture modeling functions

Σ conversion routine to full matrix

convertSigma[ {μ_, Σ_}] := Module[{type, dim, s}, type = Length[Dime ... ) & /@ Σ, 2, DiagonalMatrix[#] & /@ Σ, 3, Σ] ; s] ;

Log likelihood of single point (requires 1d data to be list)

LogLikelihood[x_, μ_, Σ_] := Module[{type, dim, s}, type = Length[Di ... [Σ], 3, Σ] ; Log[PDF[MultinormalDistribution[μ, s], x]] ] ;

Mixture likelihood of single point

MixtureLikelihoodNative[x_, {μ_, Σ_, priors_}] := Module[{type, dim, s}, typ ... ors[[i]] PDF[MultinormalDistribution[μ[[i]], s[[i]]], x], {i, 1, Length[μ]}] ] ;

Mixture log likelihood of data set

LogLikelihoodMixtureDataSet[x_, {μ_, Σ_, priors_}] := Module[{type, dim, s, ... i]]], Table[v[[j]], {j, 1, dim}]], {i, 1, Length[μ]}]] ; Plus @@ (f /@ x) ] ;

P(ω_i | θ, x_j)

ClassProbabilitiesDataSet[x_, {μ_, Σ_, priors_}] := Module[{type, dim, s, f, ... dim}]], {i, 1, Length[μ]}] ;  ((#/(Plus @@ #)) & /@ (f /@ x)) // Chop] ;

One step EM helper function

diff[x_, μ_] := (# - μ) & /@ x ;

One step of Expectation-Maximization algorithm (scalar, full covariances only)

EMOneStep[x_, {μ_, Σ_, priors_}] := Module[{ptbl, ptmp, μnew, Σnew ... amp; /@ ((ptbl tmp2) // Transpose))/ptmp] ;  {μnew, Σnew, priorsnew} ] ;

Full EM algorithm

EM[x_, {μ_, Σ_, priors_}, tst_] := NestWhileList[EMOneStep[x, #] &, {μ, Σ, priors}, (conv[#1, #2] >tst) &, 2] ;

Mixture likelihood over a range of inputs

MixtureLikelihoodMatrix[range_, n_, {μ_, Σ_, priors_}] := Module[{type, dim, ... {x, xmin, xmax + dx/2, dx}, {y, ymin, ymax + dy/2, dy}] ; Map[f, xd, {2}] // Transpose] ;

Random numbers from Gaussian mixture model

Uniform random number in range [0, 1]

RowBox[{RowBox[{RR, , :=, , RowBox[{Random, [, RowBox[{Real, ,, , RowBox[{{, RowBox[{0., ,, , 1.}], }}]}], ]}]}], ;}]

Random number generation helper functions

TransformProbabilityVector[prob_] := Module[{s = 0, len}, len = Length[prob] - 1 ; Join[Table[s = s + prob[[i]], {i, 1, len}], {1}] ] ;

PickOne[ priors_] := Module[{r}, r = RR ; (Position[(r <= #) & /@ TransformProbabilityVector[priors], True] // Flatten) // First] ;

Random number generation from Gaussian mixture model

RandomMixtureModelOne[{μ_, Σ_, priors_}] := Module[{type, dim, s, class}, ty ... ickOne[priors] ; Random[MultinormalDistribution[μ[[class]], s[[class]]]] ] ;

RandomMixtureModel[{μ_, Σ_, priors_}, n_] := Table[RandomMixtureModelOne[{μ, Σ, priors}], {n}] ;

Two-dimensional Gaussian mixture model visualization

Mixture means without overlay in 2d

RowBox[{RowBox[{PlotMeans[{μ_, Σ_, priors_}, pr_, color_:Black, opt___], :=,  ... x}, {ymin, ymax}}, FrameTrue, AspectRatioAutomatic, opt]}]}], , ]}]}], ;}]

Mixture model contour plot (with shading)

PlotContourPDF[{μ_, Σ_, priors_}, pr_, n_:50, opt___] := Module[{xmin, xmax, ... Σ, priors}], MeshRange-> {{xmin, xmax}, {ymin, ymax}}, PlotRange->All, opt] ] ;

Mixture model contour plot (no shading)

PlotContourPDFOverlay[{μ_, Σ_, priors_}, pr_, n_:50, opt___] := Module[{xmin ... hading->False, MeshRange-> {{xmin, xmax}, {ymin, ymax}}, PlotRange->All, opt] ] ;

Mixture model 3d surface plot

PlotPDF[{μ_, Σ_, priors_}, pr_, n_:50, opt___] := Module[{xmin, xmax, ymin, ... eshRange-> {{xmin, xmax}, {ymin, ymax}}, PlotRange->All, MeshFalse, opt] ] ;

Log-likelihood trajectory

PlotLogLikelihood[em_, xd_, opt___] := Module[{logp, traj, g1, g2}, logp = Log ... #62371;Show[ g1, g2, yS, AspectRatio ->1/GoldenRatio, GridLines->Automatic, opt] ] ;

Log-likelihood trajectory (varying number of classes)

PlotLogLikelihood2[em_,   xd_, opt___] := Module[{logp, traj, g1, g2, nc}, & ... #62371;Show[ g1, g2, yS, AspectRatio ->1/GoldenRatio, GridLines->Automatic, opt] ] ;

Plot multiple results for mixture model

EMAllPlot[{μ_, Σ_, priors_}, gin_, pr_: {{0, 1}, {0, 1}}] := Module[{g1, g2, ... ;False] ; Show[GraphicsArray[{{g1, gin}, {g4, g5}}], yS, ImageSize400] ] ;

Created by Mathematica  (September 8, 2003)