Perform affinity propagation clustering.

@param db Database
@param relation Relation
@return Clustering result
341053
public Clustering<MedoidModel> run(Database db, Relation<O> relation) {
    ArrayDBIDs ids = DBIDUtil.ensureArray(relation.getDBIDs());
    final int size = ids.size();

    int[] assignment = new int[size];
    double[][] s = initialization.getSimilarityMatrix(db, relation, ids);
    double[][] r = new double[size][size];
    double[][] a = new double[size][size];

    IndefiniteProgress prog = LOG.isVerbose() ? new IndefiniteProgress("Affinity Propagation Iteration", LOG) : null;
    MutableProgress aprog = LOG.isVerbose() ? new MutableProgress("Stable assignments", size + 1, LOG) : null;

    int inactive = 0;
    for(int iteration = 0; iteration < maxiter && inactive < convergence; iteration++) {
      // Update responsibility matrix:
      for(int i = 0; i < size; i++) {
        double[] ai = a[i], ri = r[i], si = s[i];
        // Find the two largest values (as initially maxk == i)
        double max1 = Double.NEGATIVE_INFINITY, max2 = Double.NEGATIVE_INFINITY;
        int maxk = -1;
        for(int k = 0; k < size; k++) {
          double val = ai[k] + si[k];
          if(val > max1) {
            max2 = max1;
            max1 = val;
            maxk = k;
          }
          else if(val > max2) {
            max2 = val;
          }
        }
        // With the maximum value known, update r:
        for(int k = 0; k < size; k++) {
          double val = si[k] - ((k != maxk) ? max1 : max2);
          ri[k] = ri[k] * lambda + val * (1. - lambda);
        }
      }
      // Update availability matrix
      for(int k = 0; k < size; k++) {
        // Compute sum of max(0, r_ik) for all i.
        // For r_kk, don't apply the max.
        double colposum = 0.;
        for(int i = 0; i < size; i++) {
          if(i == k || r[i][k] > 0.) {
            colposum += r[i][k];
          }
        }
        for(int i = 0; i < size; i++) {
          double val = colposum;
          // Adjust column sum by the one extra term.
          if(i == k || r[i][k] > 0.) {
            val -= r[i][k];
          }
          if(i != k && val > 0.) { // min
            val = 0.;
          }
          a[i][k] = a[i][k] * lambda + val * (1 - lambda);
        }
      }
      int changed = 0;
      for(int i = 0; i < size; i++) {
        double[] ai = a[i], ri = r[i];
        double max = Double.NEGATIVE_INFINITY;
        int maxj = -1;
        for(int j = 0; j < size; j++) {
          double v = ai[j] + ri[j];
          if(v > max || (i == j && v >= max)) {
            max = v;
            maxj = j;
          }
        }
        if(assignment[i] != maxj) {
          changed += 1;
          assignment[i] = maxj;
        }
      }
      inactive = (changed > 0) ? 0 : (inactive + 1);
      LOG.incrementProcessed(prog);
      if(aprog != null) {
        aprog.setProcessed(size - changed, LOG);
      }
    }
    if(aprog != null) {
      aprog.setProcessed(aprog.getTotal(), LOG);
    }
    LOG.setCompleted(prog);
    // Cluster map, by lead object
    Int2ObjectOpenHashMap<ModifiableDBIDs> map = new Int2ObjectOpenHashMap<>();
    DBIDArrayIter i1 = ids.iter();
    for(int i = 0; i1.valid(); i1.advance(), i++) {
      int c = assignment[i];
      // Add to cluster members:
      ModifiableDBIDs cids = map.get(c);
      if(cids == null) {
        cids = DBIDUtil.newArray();
        map.put(c, cids);
      }
      cids.add(i1);
    }
    // If we stopped early, the cluster lead might be in a different cluster.
    for(ObjectIterator<Int2ObjectOpenHashMap.Entry<ModifiableDBIDs>> iter = map.int2ObjectEntrySet().fastIterator(); iter.hasNext();) {
      Int2ObjectOpenHashMap.Entry<ModifiableDBIDs> entry = iter.next();
      final int key = entry.getIntKey();
      int targetkey = key;
      ModifiableDBIDs tids = null;
      // Chase arrows:
      while(tids == null && assignment[targetkey] != targetkey) {
        targetkey = assignment[targetkey];
        tids = map.get(targetkey);
      }
      if(tids != null && targetkey != key) {
        tids.addDBIDs(entry.getValue());
        iter.remove();
      }
    }

    Clustering<MedoidModel> clustering = new Clustering<>("Affinity Propagation Clustering", "ap-clustering");
    ModifiableDBIDs noise = DBIDUtil.newArray();
    for(ObjectIterator<Int2ObjectOpenHashMap.Entry<ModifiableDBIDs>> iter = map.int2ObjectEntrySet().fastIterator(); iter.hasNext();) {
      Int2ObjectOpenHashMap.Entry<ModifiableDBIDs> entry = iter.next();
      i1.seek(entry.getIntKey());
      if(entry.getValue().size() > 1) {
        MedoidModel mod = new MedoidModel(DBIDUtil.deref(i1));
        clustering.addToplevelCluster(new Cluster<>(entry.getValue(), mod));
      }
      else {
        noise.add(i1);
      }
    }
    if(noise.size() > 0) {
      MedoidModel mod = new MedoidModel(DBIDUtil.deref(noise.iter()));
      clustering.addToplevelCluster(new Cluster<>(noise, true, mod));
    }
    return clustering;
  }
337466
public Clustering<Model> run(Database database, Relation<V> relation) {
    // current dimensionality associated with each seed
    int dim_c = RelationUtil.dimensionality(relation);

    if(dim_c < l) {
      throw new IllegalStateException("Dimensionality of data < parameter l! " + "(" + dim_c + " < " + l + ")");
    }

    // current number of seeds
    int k_c = Math.min(relation.size(), k_i * k);

    // pick k0 > k points from the database
    List<ORCLUSCluster> clusters = initialSeeds(relation, k_c);

    double beta = FastMath.exp(-FastMath.log(dim_c / (double) l) * FastMath.log(1 / alpha) / FastMath.log(k_c / (double) k));

    IndefiniteProgress cprogress = LOG.isVerbose() ? new IndefiniteProgress("Current number of clusters:", LOG) : null;

    while(k_c > k) {
      // find partitioning induced by the seeds of the clusters
      assign(relation, clusters);

      // determine current subspace associated with each cluster
      for(ORCLUSCluster cluster : clusters) {
        if(cluster.objectIDs.size() > 0) {
          cluster.basis = findBasis(relation, cluster, dim_c);
        }
      }

      // reduce number of seeds and dimensionality associated with
      // each seed
      k_c = (int) Math.max(k, k_c * alpha);
      dim_c = (int) Math.max(l, dim_c * beta);
      merge(relation, clusters, k_c, dim_c, cprogress);
      if(cprogress != null) {
        cprogress.setProcessed(clusters.size(), LOG);
      }
    }
    assign(relation, clusters);

    LOG.setCompleted(cprogress);

    // get the result
    Clustering<Model> r = new Clustering<>("ORCLUS clustering", "orclus-clustering");
    for(ORCLUSCluster c : clusters) {
      r.addToplevelCluster(new Cluster<Model>(c.objectIDs, ClusterModel.CLUSTER));
    }
    return r;
  }
425016
public Result run(Database database, Relation<O> relation) {
    // Get a distance and knn query for the Euclidean distance
    // Hardcoded, only use this if you only allow the eucliden distance
    DistanceQuery<O> distQuery = database.getDistanceQuery(relation, EuclideanDistanceFunction.STATIC);
    KNNQuery<O> knnQuery = database.getKNNQuery(distQuery, 10);

    for(DBIDIter iditer = relation.iterDBIDs(); iditer.valid(); iditer.advance()) {
      // Get the actual object from the database (but discard the result)
      relation.get(iditer);
      // run a 10NN query for each point (but discard the result)
      knnQuery.getKNNForDBID(iditer, 10);
    }
    return null;
  }
865135
public OutlierResult run(Database database, Relation<O> relation) {
    final DBIDs ids = relation.getDBIDs();
    DistanceQuery<O> distFunc = database.getDistanceQuery(relation, getDistanceFunction());
    // Get k nearest neighbor and range query on the relation.
    KNNQuery<O> knnq = database.getKNNQuery(distFunc, k, DatabaseQuery.HINT_HEAVY_USE);
    RangeQuery<O> rnnQuery = database.getRangeQuery(distFunc, DatabaseQuery.HINT_HEAVY_USE);

    StepProgress stepProg = LOG.isVerbose() ? new StepProgress("DWOF", 2) : null;
    // DWOF output score storage.
    WritableDoubleDataStore dwofs = DataStoreUtil.makeDoubleStorage(ids, DataStoreFactory.HINT_DB | DataStoreFactory.HINT_HOT, 0.);
    if(stepProg != null) {
      stepProg.beginStep(1, "Initializing objects' Radii", LOG);
    }
    WritableDoubleDataStore radii = DataStoreUtil.makeDoubleStorage(ids, DataStoreFactory.HINT_TEMP | DataStoreFactory.HINT_HOT, 0.);
    // Find an initial radius for each object:
    initializeRadii(ids, knnq, distFunc, radii);
    WritableIntegerDataStore oldSizes = DataStoreUtil.makeIntegerStorage(ids, DataStoreFactory.HINT_HOT, 1);
    WritableIntegerDataStore newSizes = DataStoreUtil.makeIntegerStorage(ids, DataStoreFactory.HINT_HOT, 1);
    int countUnmerged = relation.size();
    if(stepProg != null) {
      stepProg.beginStep(2, "Clustering-Evaluating Cycles.", LOG);
    }
    IndefiniteProgress clusEvalProgress = LOG.isVerbose() ? new IndefiniteProgress("Evaluating DWOFs", LOG) : null;
    while(countUnmerged > 0) {
      LOG.incrementProcessed(clusEvalProgress);
      // Increase radii
      for(DBIDIter iter = ids.iter(); iter.valid(); iter.advance()) {
        radii.putDouble(iter, radii.doubleValue(iter) * delta);
      }
      // stores the clustering label for each object
      WritableDataStore<ModifiableDBIDs> labels = DataStoreUtil.makeStorage(ids, DataStoreFactory.HINT_TEMP, ModifiableDBIDs.class);
      // Cluster objects based on the current radius
      clusterData(ids, rnnQuery, radii, labels);
      // simple reference swap
      WritableIntegerDataStore temp = newSizes;
      newSizes = oldSizes;
      oldSizes = temp;

      // Update the cluster size count for each object.
      countUnmerged = updateSizes(ids, labels, newSizes);
      labels.destroy();
      // Update DWOF scores.
      for(DBIDIter iter = ids.iter(); iter.valid(); iter.advance()) {
        double newScore = (newSizes.intValue(iter) > 0) ? ((double) (oldSizes.intValue(iter) - 1) / (double) newSizes.intValue(iter)) : 0.0;
        dwofs.putDouble(iter, dwofs.doubleValue(iter) + newScore);
      }
    }
    LOG.setCompleted(clusEvalProgress);
    LOG.setCompleted(stepProg);
    // Build result representation.
    DoubleMinMax minmax = new DoubleMinMax();
    for(DBIDIter iter = relation.iterDBIDs(); iter.valid(); iter.advance()) {
      minmax.put(dwofs.doubleValue(iter));
    }
    OutlierScoreMeta meta = new InvertedOutlierScoreMeta(minmax.getMin(), minmax.getMax(), 0.0, Double.POSITIVE_INFINITY);
    DoubleRelation rel = new MaterializedDoubleRelation("Dynamic-Window Outlier Factors", "dwof-outlier", dwofs, ids);
    return new OutlierResult(meta, rel);
  }
190561
public double evaluateClustering(Database db, Relation<? extends NumberVector> rel, Clustering<?> c) {
    List<? extends Cluster<?>> clusters = c.getAllClusters();

    int ignorednoise = 0, withinPairs = 0;
    for(Cluster<?> cluster : clusters) {
      if((cluster.size() <= 1 || cluster.isNoise())) {
        switch(noiseHandling){
        case IGNORE_NOISE:
          ignorednoise += cluster.size();
          continue;
        case TREAT_NOISE_AS_SINGLETONS:
          continue; // No concordant distances.
        case MERGE_NOISE:
          break; // Treat like a cluster below.
        }
      }
      withinPairs += (cluster.size() * (cluster.size() - 1)) >>> 1;
      if(withinPairs < 0) {
        throw new AbortException("Integer overflow - clusters too large to compute pairwise distances.");
      }
    }
    // Materialize within-cluster distances (sorted):
    double[] withinDistances = computeWithinDistances(rel, clusters, withinPairs);
    int[] withinTies = new int[withinDistances.length];
    // Count ties within
    countTies(withinDistances, withinTies);

    long concordantPairs = 0, discordantPairs = 0, betweenPairs = 0;

    // Step two, compute discordant distances:
    for(int i = 0; i < clusters.size(); i++) {
      Cluster<?> ocluster1 = clusters.get(i);
      if((ocluster1.size() <= 1 || ocluster1.isNoise()) //
          && noiseHandling.equals(NoiseHandling.IGNORE_NOISE)) {
        continue;
      }
      for(int j = i + 1; j < clusters.size(); j++) {
        Cluster<?> ocluster2 = clusters.get(j);
        if((ocluster2.size() <= 1 || ocluster2.isNoise()) //
            && noiseHandling.equals(NoiseHandling.IGNORE_NOISE)) {
          continue;
        }
        betweenPairs += ocluster1.size() * (long) ocluster2.size();
        for(DBIDIter oit1 = ocluster1.getIDs().iter(); oit1.valid(); oit1.advance()) {
          NumberVector obj = rel.get(oit1);
          for(DBIDIter oit2 = ocluster2.getIDs().iter(); oit2.valid(); oit2.advance()) {
            double dist = distanceFunction.distance(obj, rel.get(oit2));
            int p = Arrays.binarySearch(withinDistances, dist);
            if(p >= 0) { // Tied distances:
              while(p > 0 && withinDistances[p - 1] >= dist) {
                --p;
              }
              concordantPairs += p;
              discordantPairs += withinDistances.length - p - withinTies[p];
              continue;
            }
            p = -p - 1;
            concordantPairs += p;
            discordantPairs += withinDistances.length - p;
          }
        }
      }
    }

    // Total number of pairs possible:
    final long t = ((rel.size() - ignorednoise) * (long) (rel.size() - ignorednoise - 1)) >>> 1;
    final long tt = (t * (t - 1)) >>> 1;

    double gamma = (concordantPairs - discordantPairs) / (double) (concordantPairs + discordantPairs);
    double tau = computeTau(concordantPairs, discordantPairs, tt, withinDistances.length, betweenPairs);

    // Avoid NaN when everything is in a single cluster:
    gamma = gamma > 0. ? gamma : 0.;
    tau = tau > 0. ? tau : 0.;

    if(LOG.isStatistics()) {
      LOG.statistics(new StringStatistic(key + ".noise-handling", noiseHandling.toString()));
      if(ignorednoise > 0) {
        LOG.statistics(new LongStatistic(key + ".ignored", ignorednoise));
      }
      LOG.statistics(new DoubleStatistic(key + ".gamma", gamma));
      LOG.statistics(new DoubleStatistic(key + ".tau", tau));
    }

    EvaluationResult ev = EvaluationResult.findOrCreate(db.getHierarchy(), c, "Internal Clustering Evaluation", "internal evaluation");
    MeasurementGroup g = ev.findOrCreateGroup("Concordance-based Evaluation");
    g.addMeasure("Gamma", gamma, -1., 1., 0., false);
    g.addMeasure("Tau", tau, -1., +1., 0., false);
    db.getHierarchy().resultChanged(ev);
    return gamma;
  }
1194853
public Result run(Database db, Relation<O> relation) {
    DistanceQuery<O> dq = db.getDistanceQuery(relation, getDistanceFunction());
    ArrayDBIDs ids = DBIDUtil.ensureArray(relation.getDBIDs());
    final int size = ids.size();

    if(size > 0x10000) {
      throw new AbortException("This implementation does not scale to data sets larger than " + 0x10000 + " instances (~17 GB RAM), which results in an integer overflow.");
    }
    if(Linkage.SINGLE.equals(linkage)) {
      LOG.verbose("Notice: SLINK is a much faster algorithm for single-linkage clustering!");
    }

    // Compute the initial (lower triangular) distance matrix.
    double[] scratch = new double[triangleSize(size)];
    DBIDArrayIter ix = ids.iter(), iy = ids.iter();
    // Position counter - must agree with computeOffset!
    int pos = 0;
    boolean square = Linkage.WARD.equals(linkage) && !getDistanceFunction().isSquared();
    for(int x = 0; ix.valid(); x++, ix.advance()) {
      iy.seek(0);
      for(int y = 0; y < x; y++, iy.advance()) {
        scratch[pos] = dq.distance(ix, iy);
        // Ward uses variances -- i.e. squared values
        if(square) {
          scratch[pos] *= scratch[pos];
        }
        pos++;
      }
    }

    // Initialize space for result:
    double[] height = new double[size];
    Arrays.fill(height, Double.POSITIVE_INFINITY);
    // Parent node, to track merges
    // have every object point to itself initially
    ArrayModifiableDBIDs parent = DBIDUtil.newArray(ids);
    // Active clusters, when not trivial.
    Int2ReferenceMap<ModifiableDBIDs> clusters = new Int2ReferenceOpenHashMap<>();

    // Repeat until everything merged, except the desired number of clusters:
    final int stop = size - numclusters;
    FiniteProgress prog = LOG.isVerbose() ? new FiniteProgress("Agglomerative clustering", stop, LOG) : null;
    for(int i = 0; i < stop; i++) {
      double min = Double.POSITIVE_INFINITY;
      int minx = -1, miny = -1;
      for(int x = 0; x < size; x++) {
        if(height[x] < Double.POSITIVE_INFINITY) {
          continue;
        }
        final int xbase = triangleSize(x);
        for(int y = 0; y < x; y++) {
          if(height[y] < Double.POSITIVE_INFINITY) {
            continue;
          }
          final int idx = xbase + y;
          if(scratch[idx] < min) {
            min = scratch[idx];
            minx = x;
            miny = y;
          }
        }
      }
      assert (minx >= 0 && miny >= 0);
      // Avoid allocating memory, by reusing existing iterators:
      ix.seek(minx);
      iy.seek(miny);
      // Perform merge in data structure: x -> y
      // Since y < x, prefer keeping y, dropping x.
      height[minx] = min;
      parent.set(minx, iy);
      // Merge into cluster
      ModifiableDBIDs cx = clusters.get(minx);
      ModifiableDBIDs cy = clusters.get(miny);
      int sizex = 1, sizey = 1; // cluster sizes, for averaging
      if(cy == null) {
        cy = DBIDUtil.newHashSet();
        cy.add(iy);
      }
      else {
        sizey = cy.size();
      }
      if(cx == null) {
        cy.add(ix);
      }
      else {
        sizex = cx.size();
        cy.addDBIDs(cx);
        clusters.remove(minx);
      }
      clusters.put(miny, cy);

      // Update distance matrix. Note: miny < minx

      // Implementation note: most will not need sizej, and could save the
      // hashmap lookup.
      final int xbase = triangleSize(minx), ybase = triangleSize(miny);
      // Write to (y, j), with j < y
      for(int j = 0; j < miny; j++) {
        if(height[j] < Double.POSITIVE_INFINITY) {
          continue;
        }
        final DBIDs idsj = clusters.get(j);
        final int sizej = (idsj == null) ? 1 : idsj.size();
        scratch[ybase + j] = linkage.combine(sizex, scratch[xbase + j], sizey, scratch[ybase + j], sizej, min);
      }
      // Write to (j, y), with y < j < x
      for(int j = miny + 1; j < minx; j++) {
        if(height[j] < Double.POSITIVE_INFINITY) {
          continue;
        }
        final int jbase = triangleSize(j);
        final DBIDs idsj = clusters.get(j);
        final int sizej = (idsj == null) ? 1 : idsj.size();
        scratch[jbase + miny] = linkage.combine(sizex, scratch[xbase + j], sizey, scratch[jbase + miny], sizej, min);
      }
      // Write to (j, y), with y < x < j
      for(int j = minx + 1; j < size; j++) {
        if(height[j] < Double.POSITIVE_INFINITY) {
          continue;
        }
        final DBIDs idsj = clusters.get(j);
        final int sizej = (idsj == null) ? 1 : idsj.size();
        final int jbase = triangleSize(j);
        scratch[jbase + miny] = linkage.combine(sizex, scratch[jbase + minx], sizey, scratch[jbase + miny], sizej, min);
      }
      LOG.incrementProcessed(prog);
    }
    LOG.ensureCompleted(prog);

    // Build the clustering result
    final Clustering<Model> dendrogram = new Clustering<>("Hierarchical-Clustering", "hierarchical-clustering");
    for(int x = 0; x < size; x++) {
      if(height[x] < Double.POSITIVE_INFINITY) {
        DBIDs cids = clusters.get(x);
        if(cids == null) {
          ix.seek(x);
          cids = DBIDUtil.deref(ix);
        }
        Cluster<Model> cluster = new Cluster<>("Cluster", cids);
        dendrogram.addToplevelCluster(cluster);
      }
    }

    return dendrogram;
  }
1064648
private List<Cluster<SubspaceModel>> sortClusters(Relation<V> relation, Object2ObjectMap<long[], List<ArrayModifiableDBIDs>> clustersMap) {
    final int db_dim = RelationUtil.dimensionality(relation);
    // int num = 1;
    List<Cluster<SubspaceModel>> clusters = new ArrayList<>();
    for(long[] pv : clustersMap.keySet()) {
      List<ArrayModifiableDBIDs> parallelClusters = clustersMap.get(pv);
      for(int i = 0; i < parallelClusters.size(); i++) {
        ArrayModifiableDBIDs c = parallelClusters.get(i);
        Cluster<SubspaceModel> cluster = new Cluster<>(c);
        cluster.setModel(new SubspaceModel(new Subspace(pv), Centroid.make(relation, c).getArrayRef()));
        String subspace = BitsUtil.toStringLow(cluster.getModel().getSubspace().getDimensions(), db_dim);
        cluster.setName(parallelClusters.size() > 1 ? ("Cluster_" + subspace + "_" + i) : ("Cluster_" + subspace));
        clusters.add(cluster);
      }
    }

    // sort the clusters w.r.t. lambda
    Comparator<Cluster<SubspaceModel>> comparator = new Comparator<Cluster<SubspaceModel>>() {
      @Override
      public int compare(Cluster<SubspaceModel> c1, Cluster<SubspaceModel> c2) {
        return c2.getModel().getSubspace().dimensionality() - c1.getModel().getSubspace().dimensionality();
      }
    };
    Collections.sort(clusters, comparator);
    return clusters;
  }
183650
public Clustering<CorrelationModel> run(Database database, Relation<V> relation) {
    final int dim = RelationUtil.dimensionality(relation);

    StepProgress stepprog = LOG.isVerbose() ? new StepProgress(3) : null;

    // Run Generalized DBSCAN
    LOG.beginStep(stepprog, 1, "Preprocessing local correlation dimensionalities and partitioning data");
    // FIXME: how to ensure we are running on the same relation?
    ERiCNeighborPredicate<V>.Instance npred = new ERiCNeighborPredicate<V>(settings).instantiate(database, relation);
    CorePredicate.Instance<DBIDs> cpred = new MinPtsCorePredicate(settings.minpts).instantiate(database);
    Clustering<Model> copacResult = new GeneralizedDBSCAN.Instance<>(npred, cpred, false).run();

    // extract correlation clusters
    LOG.beginStep(stepprog, 2, "Extract correlation clusters");
    List<List<Cluster<CorrelationModel>>> clusterMap = extractCorrelationClusters(copacResult, relation, dim, npred);
    if(LOG.isDebugging()) {
      StringBuilder msg = new StringBuilder("Step 2: Extract correlation clusters...");
      for(int corrDim = 0; corrDim < clusterMap.size(); corrDim++) {
        List<Cluster<CorrelationModel>> correlationClusters = clusterMap.get(corrDim);
        msg.append("\n\ncorrDim ").append(corrDim);
        for(Cluster<CorrelationModel> cluster : correlationClusters) {
          msg.append("\n  cluster ").append(cluster).append(", ids: ").append(cluster.getIDs().size());
          // .append(", level: ").append(cluster.getLevel()).append(", index:
          // ").append(cluster.getLevelIndex());
          // msg.append("\n basis " +
          // cluster.getPCA().getWeakEigenvectors().toString(" ", NF) +
          // " ids " + cluster.getIDs().size());
        }
      }
      LOG.debugFine(msg.toString());
    }
    if(LOG.isVerbose()) {
      int clusters = 0;
      for(List<Cluster<CorrelationModel>> correlationClusters : clusterMap) {
        clusters += correlationClusters.size();
      }
      LOG.verbose(clusters + " clusters extracted.");
    }

    // build hierarchy
    LOG.beginStep(stepprog, 3, "Building hierarchy");
    Clustering<CorrelationModel> clustering = new Clustering<>("ERiC clustering", "eric-clustering");
    buildHierarchy(clustering, clusterMap, npred);
    if(LOG.isDebugging()) {
      StringBuilder msg = new StringBuilder("Step 3: Build hierarchy");
      for(int corrDim = 0; corrDim < clusterMap.size(); corrDim++) {
        List<Cluster<CorrelationModel>> correlationClusters = clusterMap.get(corrDim);
        for(Cluster<CorrelationModel> cluster : correlationClusters) {
          msg.append("\n  cluster ").append(cluster).append(", ids: ").append(cluster.getIDs().size());
          // .append(", level: ").append(cluster.getLevel()).append(", index:
          // ").append(cluster.getLevelIndex());
          for(It<Cluster<CorrelationModel>> iter = clustering.getClusterHierarchy().iterParents(cluster); iter.valid(); iter.advance()) {
            msg.append("\n   parent ").append(iter.get());
          }
          for(It<Cluster<CorrelationModel>> iter = clustering.getClusterHierarchy().iterChildren(cluster); iter.valid(); iter.advance()) {
            msg.append("\n   child ").append(iter.get());
          }
        }
      }
      LOG.debugFine(msg.toString());
    }
    LOG.setCompleted(stepprog);

    for(Cluster<CorrelationModel> rc : clusterMap.get(clusterMap.size() - 1)) {
      clustering.addToplevelCluster(rc);
    }
    return clustering;
  }
391754
public OutlierResult run(Database database, Relation<N> spatial, Relation<O> relation) {
    final NeighborSetPredicate npred = getNeighborSetPredicateFactory().instantiate(database, spatial);
    DistanceQuery<O> distFunc = getNonSpatialDistanceFunction().instantiate(relation);

    WritableDoubleDataStore lrds = DataStoreUtil.makeDoubleStorage(relation.getDBIDs(), DataStoreFactory.HINT_TEMP | DataStoreFactory.HINT_HOT);
    WritableDoubleDataStore lofs = DataStoreUtil.makeDoubleStorage(relation.getDBIDs(), DataStoreFactory.HINT_STATIC);
    DoubleMinMax lofminmax = new DoubleMinMax();

    // Compute densities
    for(DBIDIter iditer = relation.iterDBIDs(); iditer.valid(); iditer.advance()) {
      DBIDs neighbors = npred.getNeighborDBIDs(iditer);
      double avg = 0;
      for(DBIDIter iter = neighbors.iter(); iter.valid(); iter.advance()) {
        avg += distFunc.distance(iditer, iter);
      }
      double lrd = 1 / (avg / neighbors.size());
      if(Double.isNaN(lrd)) {
        lrd = 0;
      }
      lrds.putDouble(iditer, lrd);
    }

    // Compute density quotients
    for(DBIDIter iditer = relation.iterDBIDs(); iditer.valid(); iditer.advance()) {
      DBIDs neighbors = npred.getNeighborDBIDs(iditer);
      double avg = 0;
      for(DBIDIter iter = neighbors.iter(); iter.valid(); iter.advance()) {
        avg += lrds.doubleValue(iter);
      }
      final double lrd = (avg / neighbors.size()) / lrds.doubleValue(iditer);
      if(!Double.isNaN(lrd)) {
        lofs.putDouble(iditer, lrd);
        lofminmax.put(lrd);
      }
      else {
        lofs.putDouble(iditer, 0.0);
      }
    }

    // Build result representation.
    DoubleRelation scoreResult = new MaterializedDoubleRelation("Spatial Outlier Factor", "sof-outlier", lofs, relation.getDBIDs());
    OutlierScoreMeta scoreMeta = new QuotientOutlierScoreMeta(lofminmax.getMin(), lofminmax.getMax(), 0.0, Double.POSITIVE_INFINITY, 1.0);
    OutlierResult or = new OutlierResult(scoreMeta, scoreResult);
    or.addChildResult(npred);
    return or;
  }
1142701
public OutlierResult run(Relation<V> relation) {
    final DBIDs ids = relation.getDBIDs();
    KNNQuery<V> knnQuery = QueryUtil.getKNNQuery(relation, getDistanceFunction(), k + 1);

    final int dim = RelationUtil.dimensionality(relation);
    if(k <= dim + 1) {
      LOG.warning("PCA is underspecified with a too low k! k should be at much larger than " + dim);
    }

    WritableDoubleDataStore cop_score = DataStoreUtil.makeDoubleStorage(ids, DataStoreFactory.HINT_HOT | DataStoreFactory.HINT_STATIC);
    WritableDataStore<double[]> cop_err_v = models ? DataStoreUtil.makeStorage(ids, DataStoreFactory.HINT_HOT | DataStoreFactory.HINT_STATIC, double[].class) : null;
    WritableIntegerDataStore cop_dim = models ? DataStoreUtil.makeIntegerStorage(ids, DataStoreFactory.HINT_HOT | DataStoreFactory.HINT_STATIC, -1) : null;
    // compute neighbors of each db object
    FiniteProgress prog = LOG.isVerbose() ? new FiniteProgress("Correlation Outlier Probabilities", relation.size(), LOG) : null;

    double[] centroid = new double[dim];
    double[] scores = new double[dim];
    ModifiableDBIDs nids = DBIDUtil.newHashSet(k + 10);
    for(DBIDIter id = ids.iter(); id.valid(); id.advance()) {
      nids.clear();
      nids.addDBIDs(knnQuery.getKNNForDBID(id, k + 1));
      nids.remove(id); // Do not use query object

      computeCentroid(centroid, relation, nids);
      PCAResult pcares = pca.processIds(nids, relation);
      double[][] tevecs = pcares.getEigenvectors();
      double[] evs = pcares.getEigenvalues();
      double[] projected = times(tevecs, minusEquals(relation.get(id).toArray(), centroid));

      if(dist == DistanceDist.CHISQUARED) {
        double sqdevs = 0;
        for(int d = 0; d < dim; d++) {
          double dev = projected[d];
          // Scale with variance and accumulate
          sqdevs += dev * dev / evs[d];
          scores[d] = 1 - ChiSquaredDistribution.cdf(sqdevs, d + 1);
        }
      }
      else {
        assert (dist == DistanceDist.GAMMA);
        double[][] dists = new double[dim][nids.size()];
        int j = 0;
        double[] srel = new double[dim];
        for(DBIDIter s = nids.iter(); s.valid() && j < nids.size(); s.advance(), j++) {
          V vec = relation.get(s);
          for(int d = 0; d < dim; d++) {
            srel[d] = vec.doubleValue(d) - centroid[d];
          }
          double sqdist = 0.0;
          for(int d = 0; d < dim; d++) {
            double serrd = transposeTimes(tevecs[d], srel);
            dists[d][j] = (sqdist += serrd * serrd / evs[d]);
          }
        }
        double sqdevs = 0;
        for(int d = 0; d < dim; d++) {
          // Scale with Stddev
          final double dev = projected[d];
          // Accumulate
          sqdevs += dev * dev / evs[d];
          // Sort, so we can trim the top 15% below.
          Arrays.sort(dists[d]);
          // Evaluate
          scores[d] = 1 - GammaChoiWetteEstimator.STATIC.estimate(dists[d], SHORTENED_ARRAY).cdf(sqdevs);
        }
      }
      // Find best score
      double min = Double.POSITIVE_INFINITY;
      int vdim = dim - 1;
      for(int d = 0; d < dim; d++) {
        double v = scores[d];
        if(v < min) {
          min = v;
          vdim = d;
        }
      }
      // Normalize the value
      final double prob = expect * (1 - min) / (expect + min);
      cop_score.putDouble(id, prob);
      if(models) {
        // Construct the error vector:
        Arrays.fill(projected, vdim + 1, dim, 0.);
        cop_err_v.put(id, timesEquals(transposeTimes(tevecs, projected), -prob));
        cop_dim.putInt(id, dim - vdim);
      }
      LOG.incrementProcessed(prog);
    }
    LOG.ensureCompleted(prog);

    // combine results.
    DoubleRelation scoreResult = new MaterializedDoubleRelation("Correlation Outlier Probabilities", COP_SCORES, cop_score, ids);
    OutlierScoreMeta scoreMeta = new ProbabilisticOutlierScore();
    OutlierResult result = new OutlierResult(scoreMeta, scoreResult);
    if(models) {
      result.addChildResult(new MaterializedRelation<>("Local Dimensionality", COP_DIM, TypeUtil.INTEGER, cop_dim, ids));
      result.addChildResult(new MaterializedRelation<>("Error vectors", COP_ERRORVEC, TypeUtil.DOUBLE_ARRAY, cop_err_v, ids));
    }
    return result;
  }
1140937
public OutlierResult run(Database database, Relation<V> relation) {
    emClustering.setSoft(true);
    Clustering<?> emresult = emClustering.run(database, relation);
    Relation<double[]> soft = null;
    for(It<Relation<double[]>> iter = emresult.getHierarchy().iterChildren(emresult).filter(Relation.class); iter.valid(); iter.advance()) {
      if(iter.get().getDataTypeInformation() == EM.SOFT_TYPE) {
        soft = iter.get();
      }
    }
    if(soft == null) {
      throw new IllegalStateException("No soft assignment generated by EM?");
    }

    double globmax = 0.0;
    WritableDoubleDataStore emo_score = DataStoreUtil.makeDoubleStorage(relation.getDBIDs(), DataStoreFactory.HINT_TEMP | DataStoreFactory.HINT_HOT);
    for(DBIDIter iditer = relation.iterDBIDs(); iditer.valid(); iditer.advance()) {
      double maxProb = Double.POSITIVE_INFINITY;
      double[] probs = soft.get(iditer);
      for(double prob : probs) {
        maxProb = Math.min(1. - prob, maxProb);
      }
      emo_score.putDouble(iditer, maxProb);
      globmax = Math.max(maxProb, globmax);
    }
    DoubleRelation scoreres = new MaterializedDoubleRelation("EM outlier scores", "em-outlier", emo_score, relation.getDBIDs());
    OutlierScoreMeta meta = new ProbabilisticOutlierScore(0.0, globmax);
    // combine results.
    OutlierResult result = new OutlierResult(meta, scoreres);
    // TODO: add a keep-EM flag?
    result.addChildResult(emresult);
    return result;
  }
268234
public OutlierResult run(Database database, Relation<O> relation) {
    DistanceQuery<O> distFunc = database.getDistanceQuery(relation, getDistanceFunction());
    RangeQuery<O> rangeQuery = database.getRangeQuery(distFunc);
    DBIDs ids = relation.getDBIDs();

    // LOCI preprocessing step
    WritableDataStore<DoubleIntArrayList> interestingDistances = DataStoreUtil.makeStorage(relation.getDBIDs(), DataStoreFactory.HINT_TEMP | DataStoreFactory.HINT_SORTED, DoubleIntArrayList.class);
    precomputeInterestingRadii(ids, rangeQuery, interestingDistances);
    // LOCI main step
    FiniteProgress progressLOCI = LOG.isVerbose() ? new FiniteProgress("LOCI scores", relation.size(), LOG) : null;
    WritableDoubleDataStore mdef_norm = DataStoreUtil.makeDoubleStorage(relation.getDBIDs(), DataStoreFactory.HINT_STATIC);
    WritableDoubleDataStore mdef_radius = DataStoreUtil.makeDoubleStorage(relation.getDBIDs(), DataStoreFactory.HINT_STATIC);
    DoubleMinMax minmax = new DoubleMinMax();

    // Shared instance, to save allocations.
    MeanVariance mv_n_r_alpha = new MeanVariance();

    for(DBIDIter iditer = ids.iter(); iditer.valid(); iditer.advance()) {
      final DoubleIntArrayList cdist = interestingDistances.get(iditer);
      final double maxdist = cdist.getDouble(cdist.size() - 1);
      final int maxneig = cdist.getInt(cdist.size() - 1);

      double maxmdefnorm = 0.0;
      double maxnormr = 0;
      if(maxneig >= nmin) {
        // Compute the largest neighborhood we will need.
        DoubleDBIDList maxneighbors = rangeQuery.getRangeForDBID(iditer, maxdist);
        // TODO: Ensure the result is sorted. This is currently implied.

        // For any critical distance, compute the normalized MDEF score.
        for(int i = 0, size = cdist.size(); i < size; i++) {
          // Only start when minimum size is fulfilled
          if(cdist.getInt(i) < nmin) {
            continue;
          }
          final double r = cdist.getDouble(i);
          final double alpha_r = alpha * r;
          // compute n(p_i, \alpha * r) from list (note: alpha_r is not cdist!)
          final int n_alphar = cdist.getInt(cdist.find(alpha_r));
          // compute \hat{n}(p_i, r, \alpha) and the corresponding \simga_{MDEF}
          mv_n_r_alpha.reset();
          for(DoubleDBIDListIter neighbor = maxneighbors.iter(); neighbor.valid(); neighbor.advance()) {
            // Stop at radius r
            if(neighbor.doubleValue() > r) {
              break;
            }
            DoubleIntArrayList cdist2 = interestingDistances.get(neighbor);
            int rn_alphar = cdist2.getInt(cdist2.find(alpha_r));
            mv_n_r_alpha.put(rn_alphar);
          }
          // We only use the average and standard deviation
          final double nhat_r_alpha = mv_n_r_alpha.getMean();
          final double sigma_nhat_r_alpha = mv_n_r_alpha.getNaiveStddev();

          // Redundant divisions by nhat_r_alpha removed.
          final double mdef = nhat_r_alpha - n_alphar;
          final double sigmamdef = sigma_nhat_r_alpha;
          final double mdefnorm = mdef / sigmamdef;

          if(mdefnorm > maxmdefnorm) {
            maxmdefnorm = mdefnorm;
            maxnormr = r;
          }
        }
      }
      else {
        // FIXME: when nmin was not fulfilled - what is the proper value then?
        maxmdefnorm = Double.POSITIVE_INFINITY;
        maxnormr = maxdist;
      }
      mdef_norm.putDouble(iditer, maxmdefnorm);
      mdef_radius.putDouble(iditer, maxnormr);
      minmax.put(maxmdefnorm);
      LOG.incrementProcessed(progressLOCI);
    }
    LOG.ensureCompleted(progressLOCI);
    DoubleRelation scoreResult = new MaterializedDoubleRelation("LOCI normalized MDEF", "loci-mdef-outlier", mdef_norm, relation.getDBIDs());
    OutlierScoreMeta scoreMeta = new QuotientOutlierScoreMeta(minmax.getMin(), minmax.getMax(), 0.0, Double.POSITIVE_INFINITY, 0.0);
    OutlierResult result = new OutlierResult(scoreMeta, scoreResult);
    result.addChildResult(new MaterializedDoubleRelation("LOCI MDEF Radius", "loci-critical-radius", mdef_radius, relation.getDBIDs()));
    return result;
  }
410705
public OutlierResult run(Relation<O> relation) {
    final DistanceQuery<O> distanceQuery = relation.getDistanceQuery(getDistanceFunction());
    final KNNQuery<O> knnQuery = relation.getKNNQuery(distanceQuery, k);

    FiniteProgress prog = LOG.isVerbose() ? new FiniteProgress("kNN distance for objects", relation.size(), LOG) : null;

    WritableDoubleDataStore knnDist = DataStoreUtil.makeDoubleStorage(relation.getDBIDs(), DataStoreFactory.HINT_HOT | DataStoreFactory.HINT_TEMP);
    WritableDBIDDataStore neighbor = DataStoreUtil.makeDBIDStorage(relation.getDBIDs(), DataStoreFactory.HINT_HOT | DataStoreFactory.HINT_TEMP);
    DBIDVar var = DBIDUtil.newVar();
    // Find nearest neighbors, and store the distances.
    for(DBIDIter it = relation.iterDBIDs(); it.valid(); it.advance()) {
      final KNNList knn = knnQuery.getKNNForDBID(it, k);
      knnDist.putDouble(it, knn.getKNNDistance());
      neighbor.put(it, knn.assignVar(knn.size() - 1, var));
      LOG.incrementProcessed(prog);
    }
    LOG.ensureCompleted(prog);

    prog = LOG.isVerbose() ? new FiniteProgress("kNN distance descriptor", relation.size(), LOG) : null;
    WritableDoubleDataStore scores = DataStoreUtil.makeDoubleStorage(relation.getDBIDs(), DataStoreFactory.HINT_DB);
    DoubleMinMax minmax = new DoubleMinMax();
    for(DBIDIter it = relation.iterDBIDs(); it.valid(); it.advance()) {
      // Distance
      double d = knnDist.doubleValue(it);
      // Distance of neighbor
      double nd = knnDist.doubleValue(neighbor.assignVar(it, var));
      double knndd = nd > 0 ? d / nd : d > 0 ? Double.POSITIVE_INFINITY : 1.;
      scores.put(it, knndd);
      minmax.put(knndd);
      LOG.incrementProcessed(prog);
    }
    LOG.ensureCompleted(prog);
    DoubleRelation scoreres = new MaterializedDoubleRelation("kNN Data Descriptor", "knndd-outlier", scores, relation.getDBIDs());
    OutlierScoreMeta meta = new BasicOutlierScoreMeta(minmax.getMin(), minmax.getMax(), 0., Double.POSITIVE_INFINITY, 1.);
    return new OutlierResult(meta, scoreres);
  }
912129
public double evaluateClustering(Database db, Relation<? extends NumberVector> rel, Clustering<?> c) {
    // FIXME: allow using a precomputed distance matrix!
    final SquaredEuclideanDistanceFunction df = SquaredEuclideanDistanceFunction.STATIC;

    List<? extends Cluster<?>> clusters = c.getAllClusters();
    double vrc = 0.;
    int ignorednoise = 0;
    if(clusters.size() > 1) {
      NumberVector[] centroids = new NumberVector[clusters.size()];
      ignorednoise = EvaluateSimplifiedSilhouette.centroids(rel, clusters, centroids, noiseOption);

      // Build global centroid and cluster count:
      final int dim = RelationUtil.dimensionality(rel);
      Centroid overallCentroid = new Centroid(dim);
      int clustercount = globalCentroid(overallCentroid, rel, clusters, centroids, noiseOption);

      // a: Distance to own centroid
      // b: Distance to overall centroid
      double a = 0, b = 0;
      Iterator<? extends Cluster<?>> ci = clusters.iterator();
      for(int i = 0; ci.hasNext(); i++) {
        Cluster<?> cluster = ci.next();
        if(cluster.size() <= 1 || cluster.isNoise()) {
          switch(noiseOption){
          case IGNORE_NOISE:
            continue; // Ignored
          case TREAT_NOISE_AS_SINGLETONS:
            // Singletons: a = 0 by definition.
            for(DBIDIter it = cluster.getIDs().iter(); it.valid(); it.advance()) {
              b += df.distance(overallCentroid, rel.get(it));
            }
            continue; // with NEXT cluster.
          case MERGE_NOISE:
            break; // Treat like a cluster below:
          }
        }
        for(DBIDIter it = cluster.getIDs().iter(); it.valid(); it.advance()) {
          NumberVector vec = rel.get(it);
          a += df.distance(centroids[i], vec);
          b += df.distance(overallCentroid, vec);
        }
      }

      vrc = ((b - a) / a) * ((rel.size() - clustercount) / (clustercount - 1.));
      // Only if {@link NoiseHandling#IGNORE_NOISE}:
      if(penalize && ignorednoise > 0) {
        vrc *= (rel.size() - ignorednoise) / (double) rel.size();
      }
    }
    if(LOG.isStatistics()) {
      LOG.statistics(new StringStatistic(key + ".vrc.noise-handling", noiseOption.toString()));
      if(ignorednoise > 0) {
        LOG.statistics(new LongStatistic(key + ".vrc.ignored", ignorednoise));
      }
      LOG.statistics(new DoubleStatistic(key + ".vrc", vrc));
    }

    EvaluationResult ev = EvaluationResult.findOrCreate(db.getHierarchy(), c, "Internal Clustering Evaluation", "internal evaluation");
    MeasurementGroup g = ev.findOrCreateGroup("Distance-based Evaluation");
    g.addMeasure("Variance Ratio Criteria", vrc, 0., 1., 0., false);
    return vrc;
  }
957243
public OutlierResult run(Relation<V> relation) {
    SimilarityQuery<V> snnInstance = similarityFunction.instantiate(relation);
    FiniteProgress progress = LOG.isVerbose() ? new FiniteProgress("Assigning Subspace Outlier Degree", relation.size(), LOG) : null;
    WritableDoubleDataStore sod_scores = DataStoreUtil.makeDoubleStorage(relation.getDBIDs(), DataStoreFactory.HINT_STATIC);
    WritableDataStore<SODModel> sod_models = models ? DataStoreUtil.makeStorage(relation.getDBIDs(), DataStoreFactory.HINT_STATIC, SODModel.class) : null;
    DoubleMinMax minmax = new DoubleMinMax();
    for(DBIDIter iter = relation.iterDBIDs(); iter.valid(); iter.advance()) {
      DBIDs neighborhood = getNearestNeighbors(relation, snnInstance, iter);

      double[] center;
      long[] weightVector = null;
      double sod = 0.;
      if(neighborhood.size() > 0) {
        center = Centroid.make(relation, neighborhood).getArrayRef();
        // Note: per-dimension variances; no covariances.
        double[] variances = computePerDimensionVariances(relation, center, neighborhood);
        double expectationOfVariance = Mean.of(variances);
        weightVector = BitsUtil.zero(variances.length);
        for(int d = 0; d < variances.length; d++) {
          if(variances[d] < alpha * expectationOfVariance) {
            BitsUtil.setI(weightVector, d);
          }
        }
        sod = subspaceOutlierDegree(relation.get(iter), center, weightVector);
      }
      else {
        center = relation.get(iter).toArray();
      }

      if(sod_models != null) {
        sod_models.put(iter, new SODModel(center, weightVector));
      }
      sod_scores.putDouble(iter, sod);
      minmax.put(sod);
      LOG.incrementProcessed(progress);
    }
    LOG.ensureCompleted(progress);
    // combine results.
    OutlierScoreMeta meta = new BasicOutlierScoreMeta(minmax.getMin(), minmax.getMax());
    OutlierResult sodResult = new OutlierResult(meta, new MaterializedDoubleRelation("Subspace Outlier Degree", "sod-outlier", sod_scores, relation.getDBIDs()));
    if(sod_models != null) {
      sodResult.addChildResult(new MaterializedRelation<>("Subspace Outlier Model", "sod-outlier", new SimpleTypeInformation<>(SODModel.class), sod_models, relation.getDBIDs()));
    }
    return sodResult;
  }
1373797
public Clustering<Model> run(Relation<?> relation) {
    HashMap<String, DBIDs> labelMap = multiple ? multipleAssignment(relation) : singleAssignment(relation);

    ModifiableDBIDs noiseids = DBIDUtil.newArray();
    Clustering<Model> result = new Clustering<>("By Label Clustering", "bylabel-clustering");
    for(Entry<String, DBIDs> entry : labelMap.entrySet()) {
      DBIDs ids = entry.getValue();
      if(ids.size() <= 1) {
        noiseids.addDBIDs(ids);
        continue;
      }
      // Build a cluster
      Cluster<Model> c = new Cluster<Model>(entry.getKey(), ids, ClusterModel.CLUSTER);
      if(noisepattern != null && noisepattern.matcher(entry.getKey()).find()) {
        c.setNoise(true);
      }
      result.addToplevelCluster(c);
    }
    // Collected noise IDs.
    if(noiseids.size() > 0) {
      Cluster<Model> c = new Cluster<Model>("Noise", noiseids, ClusterModel.CLUSTER);
      c.setNoise(true);
      result.addToplevelCluster(c);
    }
    return result;
  }
360427
@Override
  public Clustering<MeanModel> run(Database database, Relation<V> relation) {
    // Database objects to process
    final DBIDs ids = relation.getDBIDs();
    // Choose initial means
    double[][] means = initializer.chooseInitialMeans(database, relation, k, getDistanceFunction());
    // Setup cluster assignment store
    List<ModifiableDBIDs> clusters = new ArrayList<>();
    for(int i = 0; i < k; i++) {
      clusters.add(DBIDUtil.newHashSet(relation.size() / k + 2));
    }

    // Meta data storage
    final WritableDataStore<Meta> metas = initializeMeta(relation, means);
    // Perform the initial assignment
    ArrayModifiableDBIDs tids = initialAssignment(clusters, metas, ids);
    // Recompute the means after the initial assignment
    means = means(clusters, means, relation);
    // Refine the result via k-means like iterations
    means = refineResult(relation, means, clusters, metas, tids);

    // Wrap result
    Clustering<MeanModel> result = new Clustering<>("k-Means Samesize Clustering", "kmeans-samesize-clustering");
    for(int i = 0; i < clusters.size(); i++) {
      result.addToplevelCluster(new Cluster<>(clusters.get(i), new MeanModel(means[i])));
    }
    return result;
  }