Expectation Maximization

(* initialize data *)
initData:=(
Needs["MultivariateStatistics`"];
bounds=15;
nClusters=5;
npoints=10;
clusterCenters=RandomReal[{-bounds,bounds},{nClusters,2}];
data=Join@@Table[RandomReal[MultinormalDistribution[clusterCenters[[i]],IdentityMatrix[2]],{npoints}],{i,1,nClusters}]//RandomSample;
ListPlot[data]
)
(* initialize labels, centers *)
initFit:=(
labels=Join@@(Table[ConstantArray[i,npoints],{i,1,nClusters}]);
centers=RandomReal[{-bounds,bounds},{nClusters,2}];
)
mean[{}]:={0,0};
mean[data_]:=Mean[data];
(* E-step, take cluster center to be the mean of points in the cluster, update centers *)
eStep:=Module[{},centers=Table[mean[data[[Position[labels,i]//Flatten]]],{i,1,nClusters}]];
(* M-step, assign each point to closest cluster center, update labels *)
mStep:=Module[{},
distances=Array[Norm[data[[#1]]-centers[[#2]]]&,{Length[data],Length[centers]}];
labels=Table[Ordering[distances[[i]],1],{i,1,Length[data]}]//Flatten
];
initData;
initFit;

(* Visualization routine *)
visualize:=(
centerLabels={Black,Table[Text[StyleForm[i,FontSize->18,FontWeight->"Heavy"],centers[[i]]],{i,1,nClusters}]};
colors=ColorData[1,"ColorList"];
dataL=Table[{colors[[i]]}~Join~(Text[i,#]&/@(data[[Position[labels,i]//Flatten]])),{i,1,nClusters}];
Graphics[dataL~Join~centerLabels//Flatten]
);
Manipulate[Column[{
Row[{Button["E-step",eStep],Button["M-step",mStep],Button["Reset Fit",initFit],Button["Reset Data",initData]}],Dynamic[Show[visualize,PlotRange->{1.3{-bounds,bounds},{-bounds,bounds}},ImageSize->{600,600}]]}],SaveDefinitions->True]