Search in sources :

Example 81 with Plot

use of ij.gui.Plot in project GDSC-SMLM by aherbert.

the class PcPalmMolecules method calculateAveragePrecision.

/**
 * Calculate the average precision by fitting a skewed Gaussian to the histogram of the precision
 * distribution.
 *
 * <p>A simple mean and SD of the histogram is computed. If the mean of the Skewed Gaussian does
 * not fit within 3 SDs of the simple mean then the simple mean is returned.
 *
 * @param molecules the molecules
 * @param title the plot title (null if no plot should be displayed)
 * @param histogramBins the histogram bins
 * @param logFitParameters Record the fit parameters to the ImageJ log
 * @param removeOutliers The distribution is created using all values within 1.5x the
 *        inter-quartile range (IQR) of the data
 * @return The average precision
 */
public double calculateAveragePrecision(List<Molecule> molecules, String title, int histogramBins, boolean logFitParameters, boolean removeOutliers) {
    // Plot histogram of the precision
    final float[] data = new float[molecules.size()];
    final DescriptiveStatistics stats = new DescriptiveStatistics();
    double yMin = Double.NEGATIVE_INFINITY;
    double yMax = 0;
    for (int i = 0; i < data.length; i++) {
        data[i] = (float) molecules.get(i).precision;
        stats.addValue(data[i]);
    }
    // Set the min and max y-values using 1.5 x IQR
    if (removeOutliers) {
        final double lower = stats.getPercentile(25);
        final double upper = stats.getPercentile(75);
        if (Double.isNaN(lower) || Double.isNaN(upper)) {
            if (logFitParameters) {
                ImageJUtils.log("Error computing IQR: %f - %f", lower, upper);
            }
        } else {
            final double iqr = upper - lower;
            yMin = Math.max(lower - iqr, stats.getMin());
            yMax = Math.min(upper + iqr, stats.getMax());
            if (logFitParameters) {
                ImageJUtils.log("  Data range: %f - %f. Plotting 1.5x IQR: %f - %f", stats.getMin(), stats.getMax(), yMin, yMax);
            }
        }
    }
    if (yMin == Double.NEGATIVE_INFINITY) {
        yMin = stats.getMin();
        yMax = stats.getMax();
        if (logFitParameters) {
            ImageJUtils.log("  Data range: %f - %f", yMin, yMax);
        }
    }
    int bins;
    if (histogramBins <= 0) {
        bins = (int) Math.ceil((stats.getMax() - stats.getMin()) / HistogramPlot.getBinWidthScottsRule(stats.getStandardDeviation(), (int) stats.getN()));
    } else {
        bins = histogramBins;
    }
    final float[][] hist = HistogramPlot.calcHistogram(data, yMin, yMax, bins);
    Plot plot = null;
    if (title != null) {
        plot = new Plot(title, "Precision", "Frequency");
        final float[] xValues = hist[0];
        final float[] yValues = hist[1];
        plot.addPoints(xValues, yValues, Plot.BAR);
        ImageJUtils.display(title, plot);
    }
    // Extract non-zero data
    float[] x = Arrays.copyOf(hist[0], hist[0].length);
    float[] y = Arrays.copyOf(hist[1], hist[1].length);
    int count = 0;
    for (int i = 0; i < y.length; i++) {
        if (y[i] > 0) {
            x[count] = x[i];
            y[count] = y[i];
            count++;
        }
    }
    x = Arrays.copyOf(x, count);
    y = Arrays.copyOf(y, count);
    // Sense check to fitted data. Get mean and SD of histogram
    final double[] stats2 = HistogramPlot.getHistogramStatistics(x, y);
    double mean = stats2[0];
    if (logFitParameters) {
        log("  Initial Statistics: %f +/- %f", stats2[0], stats2[1]);
    }
    // Standard Gaussian fit
    final double[] parameters = fitGaussian(x, y);
    if (parameters == null) {
        log("  Failed to fit initial Gaussian");
        return mean;
    }
    double newMean = parameters[1];
    double error = Math.abs(stats2[0] - newMean) / stats2[1];
    if (error > 3) {
        log("  Failed to fit Gaussian: %f standard deviations from histogram mean", error);
        return mean;
    }
    if (newMean < yMin || newMean > yMax) {
        log("  Failed to fit Gaussian: %f outside data range %f - %f", newMean, yMin, yMax);
        return mean;
    }
    mean = newMean;
    if (logFitParameters) {
        log("  Initial Gaussian: %f @ %f +/- %f", parameters[0], parameters[1], parameters[2]);
    }
    final double[] initialSolution = new double[] { parameters[0], parameters[1], parameters[2], -1 };
    // Fit to a skewed Gaussian (or appropriate function)
    final double[] skewParameters = fitSkewGaussian(x, y, initialSolution);
    if (skewParameters == null) {
        log("  Failed to fit Skewed Gaussian");
        return mean;
    }
    final SkewNormalFunction sn = new SkewNormalFunction(skewParameters);
    if (logFitParameters) {
        log("  Skewed Gaussian: %f @ %f +/- %f (a = %f) => %f +/- %f", skewParameters[0], skewParameters[1], skewParameters[2], skewParameters[3], sn.getMean(), Math.sqrt(sn.getVariance()));
    }
    newMean = sn.getMean();
    error = Math.abs(stats2[0] - newMean) / stats2[1];
    if (error > 3) {
        log("  Failed to fit Skewed Gaussian: %f standard deviations from histogram mean", error);
        return mean;
    }
    if (newMean < yMin || newMean > yMax) {
        log("  Failed to fit Skewed Gaussian: %f outside data range %f - %f", newMean, yMin, yMax);
        return mean;
    }
    // Use original histogram x-axis to maintain all the bins
    if (plot != null) {
        plot.setColor(Color.red);
        addToPlot(plot, hist[0], skewParameters, Plot.LINE);
        plot.setColor(Color.black);
        ImageJUtils.display(title, plot);
    }
    // Return the average precision from the fitted curve
    return newMean;
}
Also used : DescriptiveStatistics(org.apache.commons.math3.stat.descriptive.DescriptiveStatistics) Plot(ij.gui.Plot) HistogramPlot(uk.ac.sussex.gdsc.core.ij.HistogramPlot) SkewNormalFunction(uk.ac.sussex.gdsc.smlm.function.SkewNormalFunction) WeightedObservedPoint(org.apache.commons.math3.fitting.WeightedObservedPoint) ClusterPoint(uk.ac.sussex.gdsc.core.clustering.ClusterPoint)

Example 82 with Plot

use of ij.gui.Plot in project GDSC-SMLM by aherbert.

the class DoubletAnalysis method plotJaccard.

private void plotJaccard(ResidualsScore residualsScore, ResidualsScore residualsScoreReference) {
    final String title = TITLE + " Score";
    final Plot plot = new Plot(title, "Residuals", "Score");
    double max = getMax(0.01, residualsScore);
    if (residualsScoreReference != null) {
        max = getMax(max, residualsScore);
    }
    plot.setLimits(0, 1, 0, max * 1.05);
    addLines(plot, residualsScore, 1);
    if (residualsScoreReference != null) {
        addLines(plot, residualsScoreReference, 0.5);
    }
    plot.setColor(Color.black);
    plot.addLabel(0, 0, String.format("Residuals %s; Jaccard %s (Black); Precision %s (Blue); Recall %s (Red)", MathUtils.rounded(residualsScore.residuals[residualsScore.maxJaccardIndex]), MathUtils.rounded(residualsScore.jaccard[residualsScore.maxJaccardIndex]), MathUtils.rounded(residualsScore.precision[residualsScore.maxJaccardIndex]), MathUtils.rounded(residualsScore.recall[residualsScore.maxJaccardIndex])));
    ImageJUtils.display(title, plot, windowOrganiser);
}
Also used : Plot(ij.gui.Plot)

Example 83 with Plot

use of ij.gui.Plot in project GDSC-SMLM by aherbert.

the class DoubletAnalysis method showCumulativeHistogram.

/**
 * Show a cumulative histogram of the data.
 *
 * @param values the values
 * @param xTitle The name of plotted statistic
 */
private void showCumulativeHistogram(double[] values, String xTitle) {
    final double[][] h = MathUtils.cumulativeHistogram(values, false);
    final String title = TITLE + " " + xTitle + " Cumulative";
    final Plot plot = new Plot(title, xTitle, "Frequency");
    plot.setColor(Color.blue);
    plot.addPoints(h[0], h[1], Plot.BAR);
    ImageJUtils.display(title, plot, windowOrganiser);
}
Also used : Plot(ij.gui.Plot)

Example 84 with Plot

use of ij.gui.Plot in project GDSC-SMLM by aherbert.

the class PsfDrift method displayPlot.

private double[][] displayPlot(String title, String yLabel, double[] x, double[] y, double[] se, LoessInterpolator loess, int start, int end) {
    // Extract non NaN numbers
    double[] newX = new double[x.length];
    double[] newY = new double[x.length];
    int count = 0;
    for (int i = 0; i < x.length; i++) {
        if (!Double.isNaN(y[i])) {
            newX[count] = x[i];
            newY[count] = y[i];
            count++;
        }
    }
    newX = Arrays.copyOf(newX, count);
    newY = Arrays.copyOf(newY, count);
    title = TITLE + " " + title;
    final Plot plot = new Plot(title, "z (nm)", yLabel);
    final double[] limitsx = MathUtils.limits(x);
    double[] limitsy = new double[2];
    if (se != null) {
        if (count > 0) {
            limitsy = new double[] { newY[0] - se[0], newY[0] + se[0] };
            for (int i = 1; i < newY.length; i++) {
                limitsy[0] = MathUtils.min(limitsy[0], newY[i] - se[i]);
                limitsy[1] = MathUtils.max(limitsy[1], newY[i] + se[i]);
            }
        }
    } else if (count > 0) {
        limitsy = MathUtils.limits(newY);
    }
    final double rangex = Math.max(0.05 * (limitsx[1] - limitsx[0]), 0.1);
    final double rangey = Math.max(0.05 * (limitsy[1] - limitsy[0]), 0.1);
    plot.setLimits(limitsx[0] - rangex, limitsx[1] + rangex, limitsy[0] - rangey, limitsy[1] + rangey);
    if (loess == null) {
        addPoints(plot, Plot.LINE, newX, newY, x[start], x[end]);
    } else {
        addPoints(plot, Plot.DOT, newX, newY, x[start], x[end]);
        newY = loess.smooth(newX, newY);
        addPoints(plot, Plot.LINE, newX, newY, x[start], x[end]);
    }
    if (se != null) {
        plot.setColor(Color.magenta);
        for (int i = 0; i < x.length; i++) {
            if (!Double.isNaN(y[i])) {
                plot.drawLine(x[i], y[i] - se[i], x[i], y[i] + se[i]);
            }
        }
        // Draw the start and end lines for the valid range
        plot.setColor(Color.green);
        plot.drawLine(x[start], limitsy[0], x[start], limitsy[1]);
        plot.drawLine(x[end], limitsy[0], x[end], limitsy[1]);
    } else {
        // draw a line for the recall limit
        plot.setColor(Color.magenta);
        plot.drawLine(limitsx[0] - rangex, settings.recallLimit, limitsx[1] + rangex, settings.recallLimit);
    }
    ImageJUtils.display(title, plot, windowOrganiser);
    return new double[][] { newX, newY };
}
Also used : Plot(ij.gui.Plot)

Example 85 with Plot

use of ij.gui.Plot in project GDSC-SMLM by aherbert.

the class PsfDrift method showHwhm.

private void showHwhm() {
    // Build a list of suitable images
    final List<String> titles = createImageList(false);
    if (titles.isEmpty()) {
        IJ.error(TITLE, "No suitable PSF images");
        return;
    }
    final GenericDialog gd = new GenericDialog(TITLE);
    gd.addMessage("Approximate the volume of the PSF as a Gaussian and\n" + "compute the equivalent Gaussian width.");
    settings = Settings.load();
    gd.addChoice("PSF", titles.toArray(new String[0]), settings.title);
    gd.addCheckbox("Use_offset", settings.useOffset);
    gd.addSlider("Smoothing", 0, 0.5, settings.smoothing);
    gd.addHelp(HelpUrls.getUrl("psf-hwhm"));
    gd.showDialog();
    if (gd.wasCanceled()) {
        return;
    }
    settings.title = gd.getNextChoice();
    settings.useOffset = gd.getNextBoolean();
    settings.smoothing = gd.getNextNumber();
    settings.save();
    imp = WindowManager.getImage(settings.title);
    if (imp == null) {
        IJ.error(TITLE, "No PSF image for image: " + settings.title);
        return;
    }
    psfSettings = getPsfSettings(imp);
    if (psfSettings == null) {
        IJ.error(TITLE, "No PSF settings for image: " + settings.title);
        return;
    }
    final int size = imp.getStackSize();
    final ImagePsfModel psf = createImagePsf(1, size, 1);
    final double[] w0 = psf.getAllHwhm0();
    final double[] w1 = psf.getAllHwhm1();
    // Get current centre
    final int centre = psfSettings.getCentreImage();
    // Extract valid values (some can be NaN)
    double[] sw0 = new double[w0.length];
    double[] sw1 = new double[w1.length];
    final TDoubleArrayList s0 = new TDoubleArrayList(w0.length);
    final TDoubleArrayList s1 = new TDoubleArrayList(w0.length);
    int c0 = 0;
    int c1 = 0;
    for (int i = 0; i < w0.length; i++) {
        if (Double.isFinite(w0[i])) {
            s0.add(i + 1);
            sw0[c0++] = w0[i];
        }
        if (Double.isFinite(w1[i])) {
            s1.add(i + 1);
            sw1[c1++] = w1[i];
        }
    }
    if (c0 == 0 && c1 == 0) {
        IJ.error(TITLE, "No computed HWHM for image: " + settings.title);
        return;
    }
    double[] slice0 = s0.toArray();
    sw0 = Arrays.copyOf(sw0, c0);
    double[] slice1 = s1.toArray();
    sw1 = Arrays.copyOf(sw1, c1);
    // Smooth
    if (settings.smoothing > 0) {
        final LoessInterpolator loess = new LoessInterpolator(settings.smoothing, 1);
        sw0 = loess.smooth(slice0, sw0);
        sw1 = loess.smooth(slice1, sw1);
    }
    final TDoubleArrayList minWx = new TDoubleArrayList();
    final TDoubleArrayList minWy = new TDoubleArrayList();
    for (int i = 0; i < w0.length; i++) {
        double weight = 0;
        if (Double.isFinite(w0[i])) {
            if (Double.isFinite(w1[i])) {
                weight = w0[i] * w1[i];
            } else {
                weight = w0[i] * w0[i];
            }
        } else if (Double.isFinite(w1[i])) {
            weight = w1[i] * w1[i];
        }
        if (weight != 0) {
            minWx.add(i + 1);
            minWy.add(Math.sqrt(weight));
        }
    }
    // Smooth the combined line
    final double[] cx = minWx.toArray();
    double[] cy = minWy.toArray();
    if (settings.smoothing > 0) {
        final LoessInterpolator loess = new LoessInterpolator(settings.smoothing, 1);
        cy = loess.smooth(cx, cy);
    }
    final int newCentre = SimpleArrayUtils.findMinIndex(cy);
    // Convert to FWHM
    final double fwhm = psfSettings.getFwhm();
    // Widths are in pixels
    final String title = TITLE + " HWHM";
    final Plot plot = new Plot(title, "Slice", "HWHM (px)");
    double[] limits = MathUtils.limits(sw0);
    limits = MathUtils.limits(limits, sw1);
    final double maxY = limits[1] * 1.05;
    plot.setLimits(1, size, 0, maxY);
    plot.setColor(Color.red);
    plot.addPoints(slice0, sw0, Plot.LINE);
    plot.setColor(Color.blue);
    plot.addPoints(slice1, sw1, Plot.LINE);
    plot.setColor(Color.magenta);
    plot.addPoints(cx, cy, Plot.LINE);
    plot.setColor(Color.black);
    plot.addLabel(0, 0, "X=red; Y=blue, Combined=Magenta");
    final PlotWindow pw = ImageJUtils.display(title, plot);
    // Show a non-blocking dialog to allow the centre to be updated ...
    // Add a label and dynamically update when the centre is moved.
    final NonBlockingExtendedGenericDialog gd2 = new NonBlockingExtendedGenericDialog(TITLE);
    final double scale = psfSettings.getPixelSize();
    // @formatter:off
    ImageJUtils.addMessage(gd2, "Update the PSF information?\n \n" + "Current z-centre = %d, FHWM = %s px (%s nm)\n", centre, MathUtils.rounded(fwhm), MathUtils.rounded(fwhm * scale));
    // @formatter:on
    gd2.addSlider("z-centre", cx[0], cx[cx.length - 1], newCentre);
    final TextField tf = gd2.getLastTextField();
    gd2.addMessage("");
    gd2.addAndGetButton("Reset", event -> tf.setText(Integer.toString(newCentre)));
    final Label label = gd2.getLastLabel();
    gd2.addCheckbox("Update_centre", settings.updateCentre);
    gd2.addCheckbox("Update_HWHM", settings.updateHwhm);
    gd2.enableYesNoCancel();
    gd2.hideCancelButton();
    final UpdateDialogListener dl = new UpdateDialogListener(cx, cy, maxY, newCentre, scale, pw, label);
    gd2.addDialogListener(dl);
    gd.addHelp(HelpUrls.getUrl("psf-hwhm"));
    gd2.showDialog();
    if (gd2.wasOKed() && (settings.updateCentre || settings.updateHwhm)) {
        final ImagePSF.Builder b = psfSettings.toBuilder();
        if (settings.updateCentre) {
            b.setCentreImage(dl.centre);
        }
        if (settings.updateHwhm) {
            b.setFwhm(dl.getFwhm());
        }
        imp.setProperty("Info", ImagePsfHelper.toString(b));
    }
}
Also used : Plot(ij.gui.Plot) NonBlockingExtendedGenericDialog(uk.ac.sussex.gdsc.core.ij.gui.NonBlockingExtendedGenericDialog) Label(java.awt.Label) PlotWindow(ij.gui.PlotWindow) LoessInterpolator(org.apache.commons.math3.analysis.interpolation.LoessInterpolator) TDoubleArrayList(gnu.trove.list.array.TDoubleArrayList) ImagePSF(uk.ac.sussex.gdsc.smlm.data.config.PSFProtos.ImagePSF) NonBlockingExtendedGenericDialog(uk.ac.sussex.gdsc.core.ij.gui.NonBlockingExtendedGenericDialog) ExtendedGenericDialog(uk.ac.sussex.gdsc.core.ij.gui.ExtendedGenericDialog) GenericDialog(ij.gui.GenericDialog) TextField(java.awt.TextField) ImagePsfModel(uk.ac.sussex.gdsc.smlm.model.ImagePsfModel)

Aggregations

Plot (ij.gui.Plot)89 HistogramPlot (uk.ac.sussex.gdsc.core.ij.HistogramPlot)20 Point (java.awt.Point)19 PlotWindow (ij.gui.PlotWindow)17 Color (java.awt.Color)13 WindowOrganiser (uk.ac.sussex.gdsc.core.ij.plugin.WindowOrganiser)13 HistogramPlotBuilder (uk.ac.sussex.gdsc.core.ij.HistogramPlot.HistogramPlotBuilder)12 BasePoint (uk.ac.sussex.gdsc.core.match.BasePoint)12 ExtendedGenericDialog (uk.ac.sussex.gdsc.core.ij.gui.ExtendedGenericDialog)11 MemoryPeakResults (uk.ac.sussex.gdsc.smlm.results.MemoryPeakResults)11 Rectangle (java.awt.Rectangle)9 ArrayList (java.util.ArrayList)9 GenericDialog (ij.gui.GenericDialog)8 NonBlockingExtendedGenericDialog (uk.ac.sussex.gdsc.core.ij.gui.NonBlockingExtendedGenericDialog)7 LocalList (uk.ac.sussex.gdsc.core.utils.LocalList)7 Statistics (uk.ac.sussex.gdsc.core.utils.Statistics)7 StoredData (uk.ac.sussex.gdsc.core.utils.StoredData)7 StoredDataStatistics (uk.ac.sussex.gdsc.core.utils.StoredDataStatistics)7 ImagePlus (ij.ImagePlus)6 TDoubleArrayList (gnu.trove.list.array.TDoubleArrayList)5