FragmentClusterer.java

Go to the documentation of this file.

package denoptim.fragmenter;
 
import java.util.ArrayList;
import java.util.Arrays;
import java.util.HashSet;
import java.util.List;
import java.util.Set;
import java.util.logging.Level;
import java.util.logging.Logger;
 
import javax.vecmath.Point3d;
 
import org.apache.commons.math3.exception.DimensionMismatchException;
import org.apache.commons.math3.ml.distance.DistanceMeasure;
import org.apache.commons.math3.stat.descriptive.SummaryStatistics;
import org.biojava.nbio.structure.geometry.CalcPoint;
import org.biojava.nbio.structure.geometry.SuperPositionSVD;
 
import denoptim.exception.DENOPTIMException;
import denoptim.graph.Fragment;
import denoptim.programs.fragmenter.FragmenterParameters;
import denoptim.utils.MathUtils;
import denoptim.utils.Randomizer;
 
public class FragmentClusterer
{   
    private List<ClusterableFragment> data;
    
    private List<DynamicCentroidCluster> clusters =
            new ArrayList<DynamicCentroidCluster>();
 
    protected FragmenterParameters settings;
    
    private Logger logger;
    
//------------------------------------------------------------------------------
    
    public FragmentClusterer(List<ClusterableFragment> data,
            FragmenterParameters settings)
    {
        this.data = data;
        this.settings = settings;
        this.logger = null;
    }
    
//------------------------------------------------------------------------------
    
    public FragmentClusterer(List<ClusterableFragment> data,
            FragmenterParameters settings, Logger logger)
    {
        this.data = data;
        this.settings = settings;
        this.logger = logger;
    }
    
//------------------------------------------------------------------------------
 
    public void cluster()
    {
        // Start by assigning each data to its own cluster
        for (int i=0; i<data.size(); i++)
        {
            DynamicCentroidCluster cluster = new DynamicCentroidCluster(data.get(i));
            clusters.add(cluster);
        }
        
        // Iteratively try to refine clusters
        boolean hasChanged = true;
        int i=0;
        while (hasChanged && i<5)
        {
            hasChanged = mergeClusters();
            i++;
        }
        
        // Print summary
        if (logger!=null)
        {
            StringBuilder sb = new StringBuilder();
            sb.append("Final number of clusters: ").append(clusters.size());
            sb.append(" (#iter. "+i+")");
            sb.append(settings.NL);
            for (int j=0; j<clusters.size(); j++)
            {
                sb.append("  Size cluster "+j+": ");
                sb.append(clusters.get(j).getPoints().size());
                sb.append(settings.NL);
            }
            logger.log(Level.INFO, sb.toString());
        }
        
    }
    
//------------------------------------------------------------------------------
    
    private boolean mergeClusters()
    {
        boolean somethingMoved = false;
        
        SuperPositionSVD svd = new SuperPositionSVD(false);
        
        Set<DynamicCentroidCluster> toRemoveClusters = 
                new HashSet<DynamicCentroidCluster>();
        for (int i=0; i<clusters.size(); i++)
        {   
            DynamicCentroidCluster clusterI = clusters.get(i);
 
            if (toRemoveClusters.contains(clusterI))
                continue;
            
            if (logger!=null)
            {
                logger.log(Level.FINE,"Clustering around centroid "+i+"...");
            }
            
            ClusterableFragment centroidI = 
                    (ClusterableFragment) clusterI.getCentroid();
            
            // Define a distance (RMSD upon superposition) for discriminating
            // this geometry from the others.
            SummaryStatistics refRMSDStats = getRMSDStatsOfNoisyDistorsions(
                    centroidI.getPoint(),
                    settings.getSizeUnimodalPop(),
                    settings.getMaxNoiseUnimodalPop());
            double rmsdThreshold = refRMSDStats.getMean() 
                    + settings.getFactorForSDOnStatsOfUnimodalPop() 
                    * refRMSDStats.getStandardDeviation();
            
            for (int j=i+1; j<clusters.size(); j++)
            {
                DynamicCentroidCluster clusterJ = 
                        clusters.get(j);
                
                if (toRemoveClusters.contains(clusterJ))
                    continue;
                
                //TODO: consider re-aligning to test alternative mappings. This
                // because the mapping is done once against the first item in 
                // the sample, but to distinguish sample members N!=1 and M!=1
                // a different mapping (i.e., a different isomorphism) might be
                // preferable.
                // Essentially, this means "get rid of the assumption that one 
                // isomorphism is suitable to align all members of the sample.
                
                ClusterableFragment centroidJ = (ClusterableFragment) 
                        clusterJ.getCentroid();
                
                Point3d[] ptsCentroidI = ClusterableFragment.convertToPointArray(
                        centroidI.getPoint());
                Point3d[] ptsCentroidJ = ClusterableFragment.convertToPointArray(
                        centroidJ.getPoint());
                svd.superposeAndTransform(ptsCentroidI, ptsCentroidJ);
                double rmsd = CalcPoint.rmsd(ptsCentroidI, ptsCentroidJ);
                if (rmsd < rmsdThreshold)
                {
                    somethingMoved = true;
                    toRemoveClusters.add(clusterJ);
                    if (logger!=null)
                    {
                        logger.log(Level.FINEST,"Merging cluster " + j + " into "
                                + "cluster " + i + " (RMSD " 
                                + String.format("%.4f", rmsd) + "<"
                                + String.format("%.4f", rmsdThreshold) + ").");
                    }
                    for (ClusterableFragment pointJ : clusterJ.getPoints())
                    {
                        Point3d[] ptsPointJ = ClusterableFragment.convertToPointArray(
                                pointJ.getPoint());
                        svd.superposeAndTransform(ptsCentroidI, ptsPointJ);
                        pointJ.setCoordsVector(ptsPointJ);
                        clusterI.addPoint(pointJ);
                    }
                } else {
                    // J looks like a cluster distinct from I. Try to move members
                    Set<ClusterableFragment> toRemoveFromJ = 
                            new HashSet<ClusterableFragment>();
                    for (ClusterableFragment pointJ : clusterJ.getPoints())
                    {
                        
                        Point3d[] ptsPointJ = ClusterableFragment.convertToPointArray(
                                pointJ.getPoint());
                        svd.superposeAndTransform(ptsCentroidI, ptsPointJ);
                        double rmsdJ = CalcPoint.rmsd(ptsCentroidJ, ptsPointJ);
                        svd.superposeAndTransform(ptsCentroidI, ptsPointJ);
                        double rmsdI = CalcPoint.rmsd(ptsCentroidI, ptsPointJ);
                        if (rmsdI < rmsdJ)
                        {
                            somethingMoved = true;
                            pointJ.setCoordsVector(ptsPointJ);
                            clusterI.addPoint(pointJ);
                            toRemoveFromJ.add(pointJ);
                            if (logger!=null)
                            {
                                logger.log(Level.FINEST,"Moving one fragment "
                                        + "from cluster " + j + " to "
                                        + "cluster " + i + " (RMSD " 
                                        + String.format("%.4f", rmsd) + ">="
                                        + String.format("%.4f", rmsdThreshold) + ").");
                            }
                        }
                    }
                    clusterJ.removeAll(toRemoveFromJ);
                    if (clusterJ.getPoints().size()==0)
                        toRemoveClusters.add(clusterJ);
                }
            }
        }
        
        clusters.removeAll(toRemoveClusters);
        
        return somethingMoved;
    }
    
//------------------------------------------------------------------------------
   
    protected static SummaryStatistics getRMSDStatsOfNoisyDistorsions(
            double[] center, int sampleSize, double maxNoise)
    {
        // We use always the same randomizer to get reproducible values 
        // no matter what randomizer is the caller using.
        Randomizer rng = new Randomizer(1L);
        
        DistanceAsRMSD measure = new DistanceAsRMSD();
        
        SummaryStatistics stats = new SummaryStatistics();
        List<double[]> refClusterCoords = new ArrayList<double[]>();
        for (int k=0; k<sampleSize; k++)
        {
            double[] coords = Arrays.copyOf(center, center.length);
            for (int j=0; j<coords.length; j++)
            {
                coords[j] = coords[j] + maxNoise*(2*rng.nextNormalDouble()-1);
            }
            refClusterCoords.add(coords);
        }
        double[] refCentroidCoords = MathUtils.centroidOf(refClusterCoords, 
                center.length);
        
        for (double[] coords : refClusterCoords) 
        {
            double rmsd = measure.compute(coords,refCentroidCoords);
            stats.addValue(rmsd);
        }
        return stats;
    }
    
//------------------------------------------------------------------------------
 
    @SuppressWarnings("serial")
    public static class DistanceAsRMSD implements DistanceMeasure
    {
        @Override
        public double compute(double[] coordsA, double[] coordsB)
                throws DimensionMismatchException
        {
            Point3d[] ptsA = new Point3d[coordsA.length/3];
            Point3d[] ptsB = new Point3d[coordsB.length/3];
            for (int i=0; i<coordsA.length/3; i++)
            {
                int j = i*3;
                ptsA[i] = new Point3d(coordsA[j],
                        coordsA[j+1],
                        coordsA[j+2]);
                ptsB[i] = new Point3d(coordsB[j],
                        coordsB[j+1],
                        coordsB[j+2]);
            }
            SuperPositionSVD svd = new SuperPositionSVD(false);
            double rmsd = svd.getRmsd(ptsA,ptsB);
            return rmsd;
        }
        
        public double compute(Point3d[] ptsA, Point3d[] ptsB)
                throws DimensionMismatchException
        {
            SuperPositionSVD svd = new SuperPositionSVD(false);
            double rmsd = svd.getRmsd(ptsA,ptsB);
            return rmsd;
        }
    }
  
//------------------------------------------------------------------------------
 
    public List<DynamicCentroidCluster> getClusters()
    {
        return clusters; 
    }
 
//------------------------------------------------------------------------------
 
    public List<List<Fragment>> getTransformedClusters()
    {
        List<List<Fragment>> transformedClusters = new ArrayList<List<Fragment>>();
        for (int i=0; i<clusters.size(); i++)
        {   
            List<Fragment> transformedCluster = new ArrayList<Fragment>();
            DynamicCentroidCluster cluster = clusters.get(i);
            for (ClusterableFragment cf : cluster.getPoints())
            {
                transformedCluster.add(cf.getTransformedCopy());
            }
            transformedClusters.add(transformedCluster);
        }
        return transformedClusters; 
    }
    
//------------------------------------------------------------------------------
 
    public List<Fragment> getClusterCentroids()
    {
        List<Fragment> centroids = new ArrayList<Fragment>();
        for (int i=0; i<clusters.size(); i++)
        {   
            DynamicCentroidCluster cluster = clusters.get(i);
            ClusterableFragment centroid = cluster.getCentroid();
            Fragment frag = centroid.getTransformedCopy();
            centroids.add(frag);
        }
        return centroids; 
    }
    
//------------------------------------------------------------------------------
    
    public List<Fragment> getNearestToClusterCentroids()
    {
        List<Fragment> nearestToCentroid = new ArrayList<Fragment>();
        for (int i=0; i<clusters.size(); i++)
        {   
            DynamicCentroidCluster cluster = clusters.get(i);
            ClusterableFragment nearest = cluster.getNearestToCentroid(
                    new DistanceAsRMSD());
            Fragment frag = nearest.getTransformedCopy();
            nearestToCentroid.add(frag);
        }
        return nearestToCentroid; 
    }
 
//------------------------------------------------------------------------------
 
}