Search in sources :

Example 1 with BinomialFitter

use of uk.ac.sussex.gdsc.smlm.fitting.BinomialFitter in project GDSC-SMLM by aherbert.

the class PcPalmClusters method fitBinomial.

/**
 * Fit a zero-truncated Binomial to the cumulative histogram.
 *
 * @param histogramData the histogram data
 * @return the double[]
 */
private double[] fitBinomial(HistogramData histogramData) {
    // Get the mean and sum of the input histogram
    double mean;
    double sum = 0;
    count = 0;
    for (int i = 0; i < histogramData.histogram[1].length; i++) {
        count += histogramData.histogram[1][i];
        sum += histogramData.histogram[1][i] * i;
    }
    mean = sum / count;
    final String name = "Zero-truncated Binomial distribution";
    ImageJUtils.log("Mean cluster size = %s", MathUtils.rounded(mean));
    ImageJUtils.log("Fitting cumulative " + name);
    // Convert to a normalised double array for the binomial fitter
    final double[] histogram = new double[histogramData.histogram[1].length];
    for (int i = 0; i < histogramData.histogram[1].length; i++) {
        histogram[i] = histogramData.histogram[1][i] / count;
    }
    // Plot the cumulative histogram
    final String title = TITLE + " Cumulative Distribution";
    Plot plot = null;
    if (settings.showCumulativeHistogram) {
        // Create a cumulative histogram for fitting
        final double[] cumulativeHistogram = new double[histogram.length];
        sum = 0;
        for (int i = 0; i < histogram.length; i++) {
            sum += histogram[i];
            cumulativeHistogram[i] = sum;
        }
        final double[] values = SimpleArrayUtils.newArray(histogram.length, 0.0, 1.0);
        plot = new Plot(title, "N", "Cumulative Probability");
        plot.addPoints(values, cumulativeHistogram, Plot.LINE);
        plot.addPoints(values, cumulativeHistogram, Plot.CIRCLE);
        // plot.addPoints(values, cumulativeHistogram, Plot.BAR);
        plot.setLimits(0, histogram.length - 1, 0, 1.05);
        ImageJUtils.display(title, plot);
    }
    // Do fitting for different N
    double bestSs = Double.POSITIVE_INFINITY;
    double[] parameters = null;
    int worse = 0;
    int countN = histogram.length - 1;
    int min = settings.minN;
    final boolean customRange = (settings.minN > 1) || (settings.maxN > 0);
    if (min > countN) {
        min = countN;
    }
    if (settings.maxN > 0 && countN > settings.maxN) {
        countN = settings.maxN;
    }
    ImageJUtils.log("Fitting N from %d to %d%s", min, countN, (customRange) ? " (custom-range)" : "");
    // Since varying the N should be done in integer steps do this
    // for n=1,2,3,... until the SS peaks then falls off (is worse then the best
    // score several times in succession)
    final BinomialFitter bf = new BinomialFitter(ImageJPluginLoggerHelper.getLogger(getClass()));
    bf.setMaximumLikelihood(settings.maximumLikelihood);
    for (int n = min; n <= countN; n++) {
        final PointValuePair solution = bf.fitBinomial(histogram, mean, n, true);
        if (solution == null) {
            continue;
        }
        final double p = solution.getPointRef()[0];
        ImageJUtils.log("Fitted %s : N=%d, p=%s. SS=%g", name, n, MathUtils.rounded(p), solution.getValue());
        if (bestSs > solution.getValue()) {
            bestSs = solution.getValue();
            parameters = new double[] { n, p };
            worse = 0;
        } else if (bestSs < Double.POSITIVE_INFINITY && ++worse >= 3) {
            break;
        }
        if (settings.showCumulativeHistogram) {
            addToPlot(n, p, title, plot, new Color((float) n / countN, 0, 1f - (float) n / countN));
        }
    }
    // Add best it in magenta
    if (settings.showCumulativeHistogram && parameters != null) {
        addToPlot((int) parameters[0], parameters[1], title, plot, Color.magenta);
    }
    return parameters;
}
Also used : Plot(ij.gui.Plot) HistogramPlot(uk.ac.sussex.gdsc.core.ij.HistogramPlot) Color(java.awt.Color) BinomialFitter(uk.ac.sussex.gdsc.smlm.fitting.BinomialFitter) ClusterPoint(uk.ac.sussex.gdsc.core.clustering.ClusterPoint) PointValuePair(org.apache.commons.math3.optim.PointValuePair)

Aggregations

Plot (ij.gui.Plot)1 Color (java.awt.Color)1 PointValuePair (org.apache.commons.math3.optim.PointValuePair)1 ClusterPoint (uk.ac.sussex.gdsc.core.clustering.ClusterPoint)1 HistogramPlot (uk.ac.sussex.gdsc.core.ij.HistogramPlot)1 BinomialFitter (uk.ac.sussex.gdsc.smlm.fitting.BinomialFitter)1