Search in sources :

Example 6 with Plot2

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

the class FIRE method showFrcTimeEvolution.

private void showFrcTimeEvolution(String name, double fireNumber, ThresholdMethod thresholdMethod, double fourierImageScale, int imageSize) {
    IJ.showStatus("Calculating FRC time evolution curve...");
    List<PeakResult> list = results.getResults();
    int nSteps = 10;
    int maxT = list.get(list.size() - 1).getFrame();
    if (maxT == 0)
        maxT = list.size();
    int step = maxT / nSteps;
    TDoubleArrayList x = new TDoubleArrayList();
    TDoubleArrayList y = new TDoubleArrayList();
    double yMin = fireNumber;
    double yMax = fireNumber;
    MemoryPeakResults newResults = new MemoryPeakResults();
    newResults.copySettings(results);
    int i = 0;
    for (int t = step; t <= maxT - step; t += step) {
        while (i < list.size()) {
            if (list.get(i).getFrame() <= t) {
                newResults.add(list.get(i));
                i++;
            } else
                break;
        }
        x.add((double) t);
        FIRE f = this.copy();
        FireResult result = f.calculateFireNumber(fourierMethod, samplingMethod, thresholdMethod, fourierImageScale, imageSize);
        double fire = (result == null) ? 0 : result.fireNumber;
        y.add(fire);
        yMin = FastMath.min(yMin, fire);
        yMax = FastMath.max(yMax, fire);
    }
    // Add the final fire number
    x.add((double) maxT);
    y.add(fireNumber);
    double[] xValues = x.toArray();
    double[] yValues = y.toArray();
    String units = "px";
    if (results.getCalibration() != null) {
        nmPerPixel = results.getNmPerPixel();
        units = "nm";
    }
    String title = name + " FRC Time Evolution";
    Plot2 plot = new Plot2(title, "Frames", "Resolution (" + units + ")", (float[]) null, (float[]) null);
    double range = Math.max(1, yMax - yMin) * 0.05;
    plot.setLimits(xValues[0], xValues[xValues.length - 1], yMin - range, yMax + range);
    plot.setColor(Color.red);
    plot.addPoints(xValues, yValues, Plot.CONNECTED_CIRCLES);
    Utils.display(title, plot);
}
Also used : TDoubleArrayList(gnu.trove.list.array.TDoubleArrayList) MemoryPeakResults(gdsc.smlm.results.MemoryPeakResults) Plot2(ij.gui.Plot2) PeakResult(gdsc.smlm.results.PeakResult) WeightedObservedPoint(org.apache.commons.math3.fitting.WeightedObservedPoint)

Example 7 with Plot2

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

the class FilterAnalysis method showPlots.

private void showPlots() {
    if (plots.isEmpty())
        return;
    // Display the top N plots
    int[] list = new int[plots.size()];
    int i = 0;
    for (NamedPlot p : plots) {
        Plot2 plot = new Plot2(p.name, p.xAxisName, "Jaccard", p.xValues, p.yValues);
        plot.setLimits(p.xValues[0], p.xValues[p.xValues.length - 1], 0, 1);
        plot.setColor(Color.RED);
        plot.draw();
        plot.setColor(Color.BLUE);
        plot.addPoints(p.xValues, p.yValues, Plot2.CROSS);
        PlotWindow plotWindow = Utils.display(p.name, plot);
        list[i++] = plotWindow.getImagePlus().getID();
    }
    new WindowOrganiser().tileWindows(list);
}
Also used : PlotWindow(ij.gui.PlotWindow) Plot2(ij.gui.Plot2) WindowOrganiser(ij.plugin.WindowOrganiser)

Example 8 with Plot2

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

the class EMGainAnalysis method fit.

/**
	 * Fit the EM-gain distribution (Gaussian * Gamma)
	 * 
	 * @param h
	 *            The distribution
	 */
private void fit(int[] h) {
    final int[] limits = limits(h);
    final double[] x = getX(limits);
    final double[] y = getY(h, limits);
    Plot2 plot = new Plot2(TITLE, "ADU", "Frequency");
    double yMax = Maths.max(y);
    plot.setLimits(limits[0], limits[1], 0, yMax);
    plot.setColor(Color.black);
    plot.addPoints(x, y, Plot2.DOT);
    Utils.display(TITLE, plot);
    // Estimate remaining parameters. 
    // Assuming a gamma_distribution(shape,scale) then mean = shape * scale
    // scale = gain
    // shape = Photons = mean / gain
    double mean = getMean(h) - bias;
    // Note: if the bias is too high then the mean will be negative. Just move the bias.
    while (mean < 0) {
        bias -= 1;
        mean += 1;
    }
    double photons = mean / gain;
    if (simulate)
        Utils.log("Simulated bias=%d, gain=%s, noise=%s, photons=%s", (int) _bias, Utils.rounded(_gain), Utils.rounded(_noise), Utils.rounded(_photons));
    Utils.log("Estimate bias=%d, gain=%s, noise=%s, photons=%s", (int) bias, Utils.rounded(gain), Utils.rounded(noise), Utils.rounded(photons));
    final int max = (int) x[x.length - 1];
    double[] g = pdf(max, photons, gain, noise, (int) bias);
    plot.setColor(Color.blue);
    plot.addPoints(x, g, Plot2.LINE);
    Utils.display(TITLE, plot);
    // Perform a fit
    CustomPowellOptimizer o = new CustomPowellOptimizer(1e-6, 1e-16, 1e-6, 1e-16);
    double[] startPoint = new double[] { photons, gain, noise, bias };
    int maxEval = 3000;
    String[] paramNames = { "Photons", "Gain", "Noise", "Bias" };
    // Set bounds
    double[] lower = new double[] { 0, 0.5 * gain, 0, bias - noise };
    double[] upper = new double[] { 2 * photons, 2 * gain, gain, bias + noise };
    // Restart until converged.
    // TODO - Maybe fix this with a better optimiser. This needs to be tested on real data.
    PointValuePair solution = null;
    for (int iter = 0; iter < 3; iter++) {
        IJ.showStatus("Fitting histogram ... Iteration " + iter);
        try {
            // Basic Powell optimiser
            MultivariateFunction fun = getFunction(limits, y, max, maxEval);
            PointValuePair optimum = o.optimize(new MaxEval(maxEval), new ObjectiveFunction(fun), GoalType.MINIMIZE, new InitialGuess((solution == null) ? startPoint : solution.getPointRef()));
            if (solution == null || optimum.getValue() < solution.getValue()) {
                double[] point = optimum.getPointRef();
                // Check the bounds
                for (int i = 0; i < point.length; i++) {
                    if (point[i] < lower[i] || point[i] > upper[i]) {
                        throw new RuntimeException(String.format("Fit out of of estimated range: %s %f", paramNames[i], point[i]));
                    }
                }
                solution = optimum;
            }
        } catch (Exception e) {
            IJ.log("Powell error: " + e.getMessage());
            if (e instanceof TooManyEvaluationsException) {
                maxEval = (int) (maxEval * 1.5);
            }
        }
        try {
            // Bounded Powell optimiser
            MultivariateFunction fun = getFunction(limits, y, max, maxEval);
            MultivariateFunctionMappingAdapter adapter = new MultivariateFunctionMappingAdapter(fun, lower, upper);
            PointValuePair optimum = o.optimize(new MaxEval(maxEval), new ObjectiveFunction(adapter), GoalType.MINIMIZE, new InitialGuess(adapter.boundedToUnbounded((solution == null) ? startPoint : solution.getPointRef())));
            double[] point = adapter.unboundedToBounded(optimum.getPointRef());
            optimum = new PointValuePair(point, optimum.getValue());
            if (solution == null || optimum.getValue() < solution.getValue()) {
                solution = optimum;
            }
        } catch (Exception e) {
            IJ.log("Bounded Powell error: " + e.getMessage());
            if (e instanceof TooManyEvaluationsException) {
                maxEval = (int) (maxEval * 1.5);
            }
        }
    }
    IJ.showStatus("");
    IJ.showProgress(1);
    if (solution == null) {
        Utils.log("Failed to fit the distribution");
        return;
    }
    double[] point = solution.getPointRef();
    photons = point[0];
    gain = point[1];
    noise = point[2];
    bias = (int) Math.round(point[3]);
    String label = String.format("Fitted bias=%d, gain=%s, noise=%s, photons=%s", (int) bias, Utils.rounded(gain), Utils.rounded(noise), Utils.rounded(photons));
    Utils.log(label);
    if (simulate) {
        Utils.log("Relative Error bias=%s, gain=%s, noise=%s, photons=%s", Utils.rounded(relativeError(bias, _bias)), Utils.rounded(relativeError(gain, _gain)), Utils.rounded(relativeError(noise, _noise)), Utils.rounded(relativeError(photons, _photons)));
    }
    // Show the PoissonGammaGaussian approximation
    double[] f = null;
    if (showApproximation) {
        f = new double[x.length];
        PoissonGammaGaussianFunction fun = new PoissonGammaGaussianFunction(1.0 / gain, noise);
        final double expected = photons * gain;
        for (int i = 0; i < f.length; i++) {
            f[i] = fun.likelihood(x[i] - bias, expected);
        //System.out.printf("x=%d, g=%f, f=%f, error=%f\n", (int) x[i], g[i], f[i],
        //		gdsc.smlm.fitting.utils.DoubleEquality.relativeError(g[i], f[i]));
        }
        yMax = Maths.maxDefault(yMax, f);
    }
    // Replot
    g = pdf(max, photons, gain, noise, (int) bias);
    plot = new Plot2(TITLE, "ADU", "Frequency");
    plot.setLimits(limits[0], limits[1], 0, yMax * 1.05);
    plot.setColor(Color.black);
    plot.addPoints(x, y, Plot2.DOT);
    plot.setColor(Color.red);
    plot.addPoints(x, g, Plot2.LINE);
    plot.addLabel(0, 0, label);
    if (showApproximation) {
        plot.setColor(Color.blue);
        plot.addPoints(x, f, Plot2.LINE);
    }
    Utils.display(TITLE, plot);
}
Also used : MaxEval(org.apache.commons.math3.optim.MaxEval) InitialGuess(org.apache.commons.math3.optim.InitialGuess) ObjectiveFunction(org.apache.commons.math3.optim.nonlinear.scalar.ObjectiveFunction) Plot2(ij.gui.Plot2) Point(java.awt.Point) TooManyEvaluationsException(org.apache.commons.math3.exception.TooManyEvaluationsException) PointValuePair(org.apache.commons.math3.optim.PointValuePair) MultivariateFunction(org.apache.commons.math3.analysis.MultivariateFunction) MultivariateFunctionMappingAdapter(org.apache.commons.math3.optim.nonlinear.scalar.MultivariateFunctionMappingAdapter) TooManyEvaluationsException(org.apache.commons.math3.exception.TooManyEvaluationsException) CustomPowellOptimizer(org.apache.commons.math3.optim.nonlinear.scalar.noderiv.CustomPowellOptimizer) PoissonGammaGaussianFunction(gdsc.smlm.function.PoissonGammaGaussianFunction)

Example 9 with Plot2

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

the class DiffusionRateTest method plotMSD.

/**
	 * Plot the MSD.
	 *
	 * @param totalSteps
	 *            the total steps
	 * @param xValues
	 *            the x values (the time)
	 * @param yValues
	 *            the y values (MSD)
	 * @param lower
	 *            the lower bounds (mean-SD)
	 * @param upper
	 *            the upper bounds (mean+SD)
	 * @param fitted
	 *            the fitted line
	 * @param dimensions
	 *            the number of dimensions for the jumps
	 */
private void plotMSD(int totalSteps, double[] xValues, double[] yValues, double[] lower, double[] upper, PolynomialFunction fitted, int dimensions) {
    // TODO Auto-generated method stub
    String title = TITLE + " " + dimensions + "D";
    Plot2 plot = new Plot2(title, "Time (seconds)", "Mean-squared Distance (um^2)", xValues, yValues);
    double[] limits = Maths.limits(upper);
    limits = Maths.limits(limits, lower);
    plot.setLimits(0, totalSteps / settings.stepsPerSecond, limits[0], limits[1]);
    plot.setColor(Color.blue);
    plot.addPoints(xValues, lower, Plot2.LINE);
    plot.addPoints(xValues, upper, Plot2.LINE);
    if (fitted != null) {
        plot.setColor(Color.red);
        plot.addPoints(new double[] { xValues[0], xValues[xValues.length - 1] }, new double[] { fitted.value(xValues[0]), fitted.value(xValues[xValues.length - 1]) }, Plot2.LINE);
    }
    plot.setColor(Color.black);
    PlotWindow pw1 = Utils.display(title, plot);
    if (Utils.isNewWindow())
        idList[idCount++] = pw1.getImagePlus().getID();
}
Also used : PlotWindow(ij.gui.PlotWindow) Plot2(ij.gui.Plot2)

Example 10 with Plot2

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

the class DensityImage method computeRipleysPlot.

/**
	 * Compute the Ripley's L-function for user selected radii and show it on a plot.
	 * 
	 * @param results
	 */
private void computeRipleysPlot(MemoryPeakResults results) {
    ExtendedGenericDialog gd = new ExtendedGenericDialog(TITLE);
    gd.addMessage("Compute Ripley's L(r) - r plot");
    gd.addNumericField("Min_radius", minR, 2);
    gd.addNumericField("Max_radius", maxR, 2);
    gd.addNumericField("Increment", incrementR, 2);
    gd.addCheckbox("Confidence_intervals", confidenceIntervals);
    gd.showDialog();
    if (gd.wasCanceled())
        return;
    minR = gd.getNextNumber();
    maxR = gd.getNextNumber();
    incrementR = gd.getNextNumber();
    confidenceIntervals = gd.getNextBoolean();
    if (minR > maxR || incrementR < 0 || gd.invalidNumber()) {
        IJ.error(TITLE, "Invalid radius parameters");
        return;
    }
    DensityManager dm = createDensityManager(results);
    double[][] values = calculateLScores(dm);
    // 99% confidence intervals
    final int iterations = (confidenceIntervals) ? 99 : 0;
    double[] upper = null;
    double[] lower = null;
    Rectangle bounds = results.getBounds();
    // Use a uniform distribution for the coordinates
    HaltonSequenceGenerator dist = new HaltonSequenceGenerator(2);
    dist.skipTo(new Well19937c(System.currentTimeMillis() + System.identityHashCode(this)).nextInt());
    for (int i = 0; i < iterations; i++) {
        IJ.showProgress(i, iterations);
        IJ.showStatus(String.format("L-score confidence interval %d / %d", i + 1, iterations));
        // Randomise coordinates
        float[] x = new float[results.size()];
        float[] y = new float[x.length];
        for (int j = x.length; j-- > 0; ) {
            final double[] d = dist.nextVector();
            x[j] = (float) (d[0] * bounds.width);
            y[j] = (float) (d[1] * bounds.height);
        }
        double[][] values2 = calculateLScores(new DensityManager(x, y, bounds));
        if (upper == null) {
            upper = values2[1];
            lower = new double[upper.length];
            System.arraycopy(upper, 0, lower, 0, upper.length);
        } else {
            for (int m = upper.length; m-- > 0; ) {
                if (upper[m] < values2[1][m])
                    upper[m] = values2[1][m];
                if (lower[m] > values2[1][m])
                    lower[m] = values2[1][m];
            }
        }
    }
    String title = results.getName() + " Ripley's (L(r) - r) / r";
    Plot2 plot = new Plot2(title, "Radius", "(L(r) - r) / r", values[0], values[1]);
    // Get the limits
    double yMin = min(0, values[1]);
    double yMax = max(0, values[1]);
    if (iterations > 0) {
        yMin = min(yMin, lower);
        yMax = max(yMax, upper);
    }
    plot.setLimits(0, values[0][values[0].length - 1], yMin, yMax);
    if (iterations > 0) {
        plot.setColor(Color.BLUE);
        plot.addPoints(values[0], upper, 1);
        plot.setColor(Color.RED);
        plot.addPoints(values[0], lower, 1);
        plot.setColor(Color.BLACK);
    }
    Utils.display(title, plot);
}
Also used : Rectangle(java.awt.Rectangle) ExtendedGenericDialog(ij.gui.ExtendedGenericDialog) DensityManager(gdsc.core.clustering.DensityManager) Plot2(ij.gui.Plot2) Well19937c(org.apache.commons.math3.random.Well19937c) HaltonSequenceGenerator(org.apache.commons.math3.random.HaltonSequenceGenerator)

Aggregations

Plot2 (ij.gui.Plot2)42 PlotWindow (ij.gui.PlotWindow)17 Point (java.awt.Point)9 BasePoint (gdsc.core.match.BasePoint)8 Statistics (gdsc.core.utils.Statistics)6 StoredDataStatistics (gdsc.core.utils.StoredDataStatistics)6 MemoryPeakResults (gdsc.smlm.results.MemoryPeakResults)5 WeightedObservedPoint (org.apache.commons.math3.fitting.WeightedObservedPoint)5 ClusterPoint (gdsc.core.clustering.ClusterPoint)4 PeakResultPoint (gdsc.smlm.ij.plugins.ResultsMatchCalculator.PeakResultPoint)4 PeakResult (gdsc.smlm.results.PeakResult)4 StoredData (gdsc.core.utils.StoredData)3 WindowOrganiser (ij.plugin.WindowOrganiser)3 Rectangle (java.awt.Rectangle)3 ArrayList (java.util.ArrayList)3 LoessInterpolator (org.apache.commons.math3.analysis.interpolation.LoessInterpolator)3 TooManyEvaluationsException (org.apache.commons.math3.exception.TooManyEvaluationsException)3 DescriptiveStatistics (org.apache.commons.math3.stat.descriptive.DescriptiveStatistics)3 Cluster (gdsc.core.clustering.Cluster)2 ClusteringEngine (gdsc.core.clustering.ClusteringEngine)2