Weka中EM算法详解

  private void EM_Init (Instances inst)
     throws Exception {
     int i, j, k;

     // 由于EM算法对初始值较敏感，故选择run k means 10 times and choose best solution
     SimpleKMeans bestK = null;
     double bestSqE = Double.MAX_VALUE;
     for (i = 0; i < 10; i++) {
       SimpleKMeans sk = new SimpleKMeans();
       sk.setSeed(m_rr.nextInt());
       sk.setNumClusters(m_num_clusters);
       sk.setDisplayStdDevs(true);
       sk.buildClusterer(inst);
       //KMeans中各个cluster的平方误差
       if (sk.getSquaredError() < bestSqE) {

           bestSqE = sk.getSquaredError();
           bestK = sk;
       }
     }

     /*************** KMeans Finds the best cluster number *****************/

     // initialize with best k-means solution
     m_num_clusters = bestK.numberOfClusters();
     // 每个样本所在各个集群的概率
     m_weights = new double[inst.numInstances()][m_num_clusters];
     // 评估每个集群所对应的离散型属性的相关取值
30     m_model = new DiscreteEstimator[m_num_clusters][m_num_attribs];
     // 每个集群所对应的连续性属性数所对应的相关取值（均值，标准偏差，样本权值（进行归一化））
m_modelNormal = new double[m_num_clusters][m_num_attribs][3];
     // 每个集群所对应的先验概率
     m_priors = new double[m_num_clusters];
     // 每个集群所对应的中心点
     Instances centers = bestK.getClusterCentroids();
     // 每个集群所对应的标准差
     Instances stdD = bestK.getClusterStandardDevs();
     // ??? Returns for each cluster the frequency counts for the values of each nominal attribute
     int [][][] nominalCounts = bestK.getClusterNominalCounts();
     // 得到每个集群所对应的样本数
     int [] clusterSizes = bestK.getClusterSizes();

     for (i = 0; i < m_num_clusters; i++) {
       Instance center = centers.instance(i);
       for (j = 0; j < m_num_attribs; j++) {

       // 样本属性是离散型
       if (inst.attribute(j).isNominal())
       {
           m_model[i][j] = new DiscreteEstimator(m_theInstances.attribute(j).numValues()
                                               , true);
           for (k = 0; k < inst.attribute(j).numValues(); k++) {
               m_model[i][j].addValue(k, nominalCounts[i][j][k]);
           }
       }
       //// 样本属性是连续型
       else
       {
           double minStdD = (m_minStdDevPerAtt != null)? m_minStdDevPerAtt[j]: m_minStdDev;
           double mean = (center.isMissing(j))? inst.meanOrMode(j): center.value(j);
           m_modelNormal[i][j][0] = mean;
           double stdv = (stdD.instance(i).isMissing(j))? ((m_maxValues[j] -
                       m_minValues[j]) / (2 * m_num_clusters)): stdD.instance(i).value(j);
           if (stdv < minStdD)
           {
               stdv = inst.attributeStats(j).numericStats.stdDev;
               if (Double.isInfinite(stdv)) {
                   stdv = minStdD;
               }
               if (stdv < minStdD) {
                   stdv = minStdD;
               }
           }
           if (stdv <= 0) {
               stdv = m_minStdDev;
           }

           m_modelNormal[i][j][1] = stdv;
           m_modelNormal[i][j][2] = 1.0;
       }
      }
     }    

     for (j = 0; j < m_num_clusters; j++) {
       // 计算每个集群的先验概率
       m_priors[j] = clusterSizes[j];
     }
     Utils.normalize(m_priors);
   }