EM算法求解高斯混合模型

public class GausssianEM {
public static void main(String[] args) {
int k=4,n=10,i=1;
Gaussian[] gaussians=new Gaussian[k];
double[][] R=new double
[k];
gaussians[0]=new Gaussian(0.1,2,1);
gaussians[1]=new Gaussian(0.4,0,0.4);
gaussians[2]=new Gaussian(0.3,3,2);
gaussians[3]=new Gaussian(0.2,-1,6);
Gaussian[] newGs=gaussians;
double[] Y={1,-0.3,4.5,2,1.2,8,10,-2.7,3.1,5};
while(i<=1000){
R=Estep(Y,gaussians);
newGs=Mstep(R,Y,gaussians);
if(distance(newGs, gaussians)<0.0000001){
gaussians=newGs;break;
}
gaussians=newGs;i++;
}
System.out.println("迭代次数:"+i);
for(Gaussian gs:gaussians)
System.out.println(gs.toString());
}
public static double[][] Estep(double[] Y,Gaussian[] gaussians){
double[][] R=new double[Y.length][gaussians.length];
double sum;
for(int j=0;j<Y.length;j++){
sum=0;
for(int k=0;k<gaussians.length;k++){
R[j][k]=gaussians[k].valueAt(Y[j]);
sum+=R[j][k];
}
for(int k=0;k<gaussians.length;k++) R[j][k]=R[j][k]/sum;
}
return R;
}
public static Gaussian[] Mstep(double[][] R,double[] Y,Gaussian[] gaussians){
Gaussian[] newGs=new Gaussian[gaussians.length];
double sumRk,sumRk_Y,sumRk_YU;
for(int k=0;k<gaussians.length;k++){
sumRk=0;sumRk_Y=0;sumRk_YU=0;
for(int j=0;j<Y.length;j++){
sumRk+=R[j][k];
sumRk_Y+=R[j][k]*Y[j];
sumRk_YU+=R[j][k]*Math.pow(Y[j]-gaussians[k].avgValue, 2);
}
newGs[k]=new Gaussian(sumRk/Y.length, sumRk_Y/sumRk, sumRk_YU/sumRk);
}
return newGs;
}
public static double distance(Gaussian[] news,Gaussian[] olds){
double distance=0;
for(int i=0;i<olds.length;i++){
distance+=Math.pow(news[i].weight-olds[i].weight, 2);
distance+=Math.pow(news[i].var-olds[i].var, 2);
distance+=Math.pow(news[i].avgValue-olds[i].avgValue, 2);
}
return Math.pow(distance, 0.5);
}
}
class Gaussian{
double weight,var,avgValue;
public Gaussian(double weight, double avgValue,double var) {
this.weight = weight;this.var = var;this.avgValue = avgValue;
}
public double valueAt(double y){
double result=Math.exp(-(double)Math.pow((y-avgValue),2)/(2*var));
result=(1.0/(Math.pow(6.28*var,0.5)))*result;
return result;
}
public String toString() {
return "Gaussian [weight=" + weight + ", var=" + var + ", avgValue="
+ avgValue + "]";
}
}