Examples with Mean org.apache.commons.math3.stat.descriptive.moment.Mean used on opensource projects

Example 6 with Mean

use of org.apache.commons.math3.stat.descriptive.moment.Mean in project GDSC-SMLM by aherbert.

the class CreateData method createPhotonDistribution.

/**
	 * @return A photon distribution loaded from a file of floating-point values with the specified population mean.
	 */
private RealDistribution createPhotonDistribution() {
    if (PHOTON_DISTRIBUTION[PHOTON_CUSTOM].equals(settings.photonDistribution)) {
        // Get the distribution file
        String filename = Utils.getFilename("Photon_distribution", settings.photonDistributionFile);
        if (filename != null) {
            settings.photonDistributionFile = filename;
            try {
                InputStream is = new FileInputStream(new File(settings.photonDistributionFile));
                BufferedReader in = new BufferedReader(new UnicodeReader(is, null));
                StoredDataStatistics stats = new StoredDataStatistics();
                try {
                    String str = null;
                    double val = 0.0d;
                    while ((str = in.readLine()) != null) {
                        val = Double.parseDouble(str);
                        stats.add(val);
                    }
                } finally {
                    in.close();
                }
                if (stats.getSum() > 0) {
                    // Update the statistics to the desired mean.
                    double scale = (double) settings.photonsPerSecond / stats.getMean();
                    double[] values = stats.getValues();
                    for (int i = 0; i < values.length; i++) values[i] *= scale;
                    // TODO - Investigate the limits of this distribution. 
                    // How far above and below the input data will values be generated.
                    // Create the distribution using the recommended number of bins
                    final int binCount = stats.getN() / 10;
                    EmpiricalDistribution dist = new EmpiricalDistribution(binCount, createRandomGenerator());
                    dist.load(values);
                    return dist;
                }
            } catch (IOException e) {
            // Ignore
            } catch (NullArgumentException e) {
            // Ignore 
            } catch (NumberFormatException e) {
            // Ignore
            }
        }
        Utils.log("Failed to load custom photon distribution from file: %s. Default to fixed.", settings.photonDistributionFile);
    } else if (PHOTON_DISTRIBUTION[PHOTON_UNIFORM].equals(settings.photonDistribution)) {
        if (settings.photonsPerSecond < settings.photonsPerSecondMaximum) {
            UniformRealDistribution dist = new UniformRealDistribution(createRandomGenerator(), settings.photonsPerSecond, settings.photonsPerSecondMaximum);
            return dist;
        }
    } else if (PHOTON_DISTRIBUTION[PHOTON_GAMMA].equals(settings.photonDistribution)) {
        final double scaleParameter = settings.photonsPerSecond / settings.photonShape;
        GammaDistribution dist = new GammaDistribution(createRandomGenerator(), settings.photonShape, scaleParameter, ExponentialDistribution.DEFAULT_INVERSE_ABSOLUTE_ACCURACY);
        return dist;
    } else if (PHOTON_DISTRIBUTION[PHOTON_CORRELATED].equals(settings.photonDistribution)) {
        // No distribution required
        return null;
    }
    settings.photonDistribution = PHOTON_DISTRIBUTION[PHOTON_FIXED];
    return null;
}

Example 7 with Mean

use of org.apache.commons.math3.stat.descriptive.moment.Mean in project GDSC-SMLM by aherbert.

the class CMOSAnalysis method computeError.

private void computeError(int slice, ImageStack simulationStack) {
    String label = simulationStack.getSliceLabel(slice);
    float[] e = (float[]) simulationStack.getPixels(slice);
    float[] o = (float[]) measuredStack.getPixels(slice);
    // Get the mean error
    Statistics s = new Statistics();
    for (int i = e.length; i-- > 0; ) s.add(o[i] - e[i]);
    StringBuilder result = new StringBuilder("Error ").append(label);
    result.append(" = ").append(Utils.rounded(s.getMean()));
    result.append(" +/- ").append(Utils.rounded(s.getStandardDeviation()));
    // Do statistical tests
    double[] x = Utils.toDouble(e), y = Utils.toDouble(o);
    PearsonsCorrelation c = new PearsonsCorrelation();
    result.append(" : R=").append(Utils.rounded(c.correlation(x, y)));
    // Mann-Whitney U is valid for any distribution, e.g. variance
    MannWhitneyUTest test = new MannWhitneyUTest();
    double p = test.mannWhitneyUTest(x, y);
    result.append(" : Mann-Whitney U p=").append(Utils.rounded(p)).append(' ').append(((p < 0.05) ? "reject" : "accept"));
    if (slice != 2) {
        // T-Test is valid for approximately Normal distributions, e.g. offset and gain
        p = TestUtils.tTest(x, y);
        result.append(" : T-Test p=").append(Utils.rounded(p)).append(' ').append(((p < 0.05) ? "reject" : "accept"));
        p = TestUtils.pairedTTest(x, y);
        result.append(" : Paired T-Test p=").append(Utils.rounded(p)).append(' ').append(((p < 0.05) ? "reject" : "accept"));
    }
    Utils.log(result.toString());
}

Example 8 with Mean

use of org.apache.commons.math3.stat.descriptive.moment.Mean in project GDSC-SMLM by aherbert.

the class BaseFunctionSolverTest method canFitSingleGaussianBetter.

void canFitSingleGaussianBetter(FunctionSolver solver, boolean applyBounds, FunctionSolver solver2, boolean applyBounds2, String name, String name2, NoiseModel noiseModel) {
    double[] noise = getNoise(noiseModel);
    if (solver.isWeighted())
        solver.setWeights(getWeights(noiseModel));
    int LOOPS = 5;
    randomGenerator.setSeed(seed);
    StoredDataStatistics[] stats = new StoredDataStatistics[6];
    String[] statName = { "Signal", "X", "Y" };
    int[] betterPrecision = new int[3];
    int[] totalPrecision = new int[3];
    int[] betterAccuracy = new int[3];
    int[] totalAccuracy = new int[3];
    int i1 = 0, i2 = 0;
    for (double s : signal) {
        double[] expected = createParams(1, s, 0, 0, 1);
        double[] lower = null, upper = null;
        if (applyBounds || applyBounds2) {
            lower = createParams(0, s * 0.5, -0.2, -0.2, 0.8);
            upper = createParams(3, s * 2, 0.2, 0.2, 1.2);
        }
        if (applyBounds)
            solver.setBounds(lower, upper);
        if (applyBounds2)
            solver2.setBounds(lower, upper);
        for (int loop = LOOPS; loop-- > 0; ) {
            double[] data = drawGaussian(expected, noise, noiseModel);
            for (int i = 0; i < stats.length; i++) stats[i] = new StoredDataStatistics();
            for (double db : base) for (double dx : shift) for (double dy : shift) for (double dsx : factor) {
                double[] p = createParams(db, s, dx, dy, dsx);
                double[] fp = fitGaussian(solver, data, p, expected);
                i1 += solver.getEvaluations();
                double[] fp2 = fitGaussian(solver2, data, p, expected);
                i2 += solver2.getEvaluations();
                // Get the mean and sd (the fit precision)
                compare(fp, expected, fp2, expected, Gaussian2DFunction.SIGNAL, stats[0], stats[1]);
                compare(fp, expected, fp2, expected, Gaussian2DFunction.X_POSITION, stats[2], stats[3]);
                compare(fp, expected, fp2, expected, Gaussian2DFunction.Y_POSITION, stats[4], stats[5]);
            // Use the distance
            //stats[2].add(distance(fp, expected));
            //stats[3].add(distance(fp2, expected2));
            }
            // two sided
            double alpha = 0.05;
            for (int i = 0; i < stats.length; i += 2) {
                double u1 = stats[i].getMean();
                double u2 = stats[i + 1].getMean();
                double sd1 = stats[i].getStandardDeviation();
                double sd2 = stats[i + 1].getStandardDeviation();
                TTest tt = new TTest();
                boolean diff = tt.tTest(stats[i].getValues(), stats[i + 1].getValues(), alpha);
                int index = i / 2;
                String msg = String.format("%s vs %s : %.1f (%s) %s %f +/- %f vs %f +/- %f  (N=%d) %b", name2, name, s, noiseModel, statName[index], u2, sd2, u1, sd1, stats[i].getN(), diff);
                if (diff) {
                    // Different means. Check they are roughly the same
                    if (DoubleEquality.almostEqualRelativeOrAbsolute(u1, u2, 0.1, 0)) {
                        // Basically the same. Check which is more precise
                        if (!DoubleEquality.almostEqualRelativeOrAbsolute(sd1, sd2, 0.05, 0)) {
                            if (sd2 < sd1) {
                                betterPrecision[index]++;
                                println(msg + " P*");
                            } else
                                println(msg + " P");
                            totalPrecision[index]++;
                        }
                    } else {
                        // Check which is more accurate (closer to zero)
                        u1 = Math.abs(u1);
                        u2 = Math.abs(u2);
                        if (u2 < u1) {
                            betterAccuracy[index]++;
                            println(msg + " A*");
                        } else
                            println(msg + " A");
                        totalAccuracy[index]++;
                    }
                } else {
                    // The same means. Check that it is more precise
                    if (!DoubleEquality.almostEqualRelativeOrAbsolute(sd1, sd2, 0.05, 0)) {
                        if (sd2 < sd1) {
                            betterPrecision[index]++;
                            println(msg + " P*");
                        } else
                            println(msg + " P");
                        totalPrecision[index]++;
                    }
                }
            }
        }
    }
    int better = 0, total = 0;
    for (int index = 0; index < statName.length; index++) {
        better += betterPrecision[index] + betterAccuracy[index];
        total += totalPrecision[index] + totalAccuracy[index];
        test(name2, name, statName[index] + " P", betterPrecision[index], totalPrecision[index], printBetterDetails);
        test(name2, name, statName[index] + " A", betterAccuracy[index], totalAccuracy[index], printBetterDetails);
    }
    test(name2, name, String.format("All (eval [%d] [%d]) : ", i2, i1), better, total, true);
}

Example 9 with Mean

use of org.apache.commons.math3.stat.descriptive.moment.Mean in project GDSC-SMLM by aherbert.

the class FIRE method showDialog.

private boolean showDialog() {
    ExtendedGenericDialog gd = new ExtendedGenericDialog(TITLE);
    gd.addMessage("Compute the resolution using Fourier Ring Correlation");
    gd.addHelp(About.HELP_URL);
    boolean single = results2 == null;
    gd.addMessage("Image construction options:");
    gd.addChoice("Image_scale", SCALE_ITEMS, SCALE_ITEMS[imageScaleIndex]);
    gd.addChoice("Auto_image_size", IMAGE_SIZE_ITEMS, IMAGE_SIZE_ITEMS[imageSizeIndex]);
    if (extraOptions)
        gd.addCheckbox("Use_signal (if present)", useSignal);
    gd.addNumericField("Max_per_bin", maxPerBin, 0);
    gd.addMessage("Fourier options:");
    String[] fourierMethodNames = SettingsManager.getNames((Object[]) FRC.FourierMethod.values());
    gd.addChoice("Fourier_method", fourierMethodNames, fourierMethodNames[fourierMethodIndex]);
    String[] samplingMethodNames = SettingsManager.getNames((Object[]) FRC.SamplingMethod.values());
    gd.addChoice("Sampling_method", samplingMethodNames, samplingMethodNames[samplingMethodIndex]);
    gd.addSlider("Sampling_factor", 0.2, 4, perimeterSamplingFactor);
    gd.addMessage("FIRE options:");
    String[] thresholdMethodNames = SettingsManager.getNames((Object[]) FRC.ThresholdMethod.values());
    gd.addChoice("Threshold_method", thresholdMethodNames, thresholdMethodNames[thresholdMethodIndex]);
    gd.addCheckbox("Show_FRC_curve", showFRCCurve);
    if (single) {
        gd.addMessage("For single datasets:");
        Label l = (Label) gd.getMessage();
        gd.addNumericField("Block_size", blockSize, 0);
        gd.addCheckbox("Random_split", randomSplit);
        gd.addNumericField("Repeats", repeats, 0);
        gd.addCheckbox("Show_FRC_curve_repeats", showFRCCurveRepeats);
        gd.addCheckbox("Show_FRC_time_evolution", showFRCTimeEvolution);
        gd.addCheckbox("Spurious correlation correction", spuriousCorrelationCorrection);
        gd.addNumericField("Q-value", qValue, 3);
        gd.addNumericField("Precision_Mean", mean, 2, 6, "nm");
        gd.addNumericField("Precision_Sigma", sigma, 2, 6, "nm");
        if (extraOptions)
            gd.addNumericField("Threads", getLastNThreads(), 0);
        // Rearrange the dialog
        if (gd.getLayout() != null) {
            GridBagLayout grid = (GridBagLayout) gd.getLayout();
            int xOffset = 0, yOffset = 0;
            int lastY = -1, rowCount = 0;
            for (Component comp : gd.getComponents()) {
                // Check if this should be the second major column
                if (comp == l) {
                    xOffset += 2;
                    // Skip title row
                    yOffset = yOffset - rowCount + 1;
                }
                // Reposition the field
                GridBagConstraints c = grid.getConstraints(comp);
                if (lastY != c.gridy)
                    rowCount++;
                lastY = c.gridy;
                c.gridx = c.gridx + xOffset;
                c.gridy = c.gridy + yOffset;
                c.insets.left = c.insets.left + 10 * xOffset;
                c.insets.top = 0;
                c.insets.bottom = 0;
                grid.setConstraints(comp, c);
            }
            if (IJ.isLinux())
                gd.setBackground(new Color(238, 238, 238));
        }
    }
    gd.showDialog();
    if (gd.wasCanceled())
        return false;
    imageScaleIndex = gd.getNextChoiceIndex();
    imageSizeIndex = gd.getNextChoiceIndex();
    if (extraOptions)
        myUseSignal = useSignal = gd.getNextBoolean();
    maxPerBin = Math.abs((int) gd.getNextNumber());
    fourierMethodIndex = gd.getNextChoiceIndex();
    fourierMethod = FourierMethod.values()[fourierMethodIndex];
    samplingMethodIndex = gd.getNextChoiceIndex();
    samplingMethod = SamplingMethod.values()[samplingMethodIndex];
    perimeterSamplingFactor = gd.getNextNumber();
    thresholdMethodIndex = gd.getNextChoiceIndex();
    thresholdMethod = FRC.ThresholdMethod.values()[thresholdMethodIndex];
    showFRCCurve = gd.getNextBoolean();
    if (single) {
        blockSize = Math.max(1, (int) gd.getNextNumber());
        randomSplit = gd.getNextBoolean();
        repeats = Math.max(1, (int) gd.getNextNumber());
        showFRCCurveRepeats = gd.getNextBoolean();
        showFRCTimeEvolution = gd.getNextBoolean();
        spuriousCorrelationCorrection = gd.getNextBoolean();
        qValue = Math.abs(gd.getNextNumber());
        mean = Math.abs(gd.getNextNumber());
        sigma = Math.abs(gd.getNextNumber());
        if (extraOptions) {
            setThreads((int) gd.getNextNumber());
            lastNThreads = this.nThreads;
        }
    }
    // Check arguments
    try {
        Parameters.isAboveZero("Perimeter sampling factor", perimeterSamplingFactor);
        if (single && spuriousCorrelationCorrection) {
            Parameters.isAboveZero("Q-value", qValue);
            Parameters.isAboveZero("Precision Mean", mean);
            Parameters.isAboveZero("Precision Sigma", sigma);
            // Set these for use in FIRE computation 
            setCorrectionParameters(qValue, mean, sigma);
        }
    } catch (IllegalArgumentException e) {
        IJ.error(TITLE, e.getMessage());
        return false;
    }
    return true;
}

Example 10 with Mean

use of org.apache.commons.math3.stat.descriptive.moment.Mean in project GDSC-SMLM by aherbert.

the class FIRE method showQEstimationInputDialog.

private boolean showQEstimationInputDialog() {
    ExtendedGenericDialog gd = new ExtendedGenericDialog(TITLE);
    gd.addHelp(About.HELP_URL);
    // Build a list of all images with a region ROI
    List<String> titles = new LinkedList<String>();
    if (WindowManager.getWindowCount() > 0) {
        for (int imageID : WindowManager.getIDList()) {
            ImagePlus imp = WindowManager.getImage(imageID);
            if (imp != null && imp.getRoi() != null && imp.getRoi().isArea())
                titles.add(imp.getTitle());
        }
    }
    gd.addMessage("Estimate the blinking correction parameter Q for Fourier Ring Correlation");
    ResultsManager.addInput(gd, inputOption, InputSource.MEMORY);
    if (!titles.isEmpty())
        gd.addCheckbox((titles.size() == 1) ? "Use_ROI" : "Choose_ROI", chooseRoi);
    gd.addMessage("Image construction options:");
    //gd.addCheckbox("Use_signal (if present)", useSignal);
    gd.addChoice("Image_scale", SCALE_ITEMS, SCALE_ITEMS[imageScaleIndex]);
    gd.addChoice("Auto_image_size", IMAGE_SIZE_ITEMS, IMAGE_SIZE_ITEMS[imageSizeIndex]);
    gd.addNumericField("Block_size", blockSize, 0);
    gd.addCheckbox("Random_split", randomSplit);
    gd.addNumericField("Max_per_bin", maxPerBin, 0);
    gd.addMessage("Fourier options:");
    String[] fourierMethodNames = SettingsManager.getNames((Object[]) FRC.FourierMethod.values());
    gd.addChoice("Fourier_method", fourierMethodNames, fourierMethodNames[fourierMethodIndex]);
    String[] samplingMethodNames = SettingsManager.getNames((Object[]) FRC.SamplingMethod.values());
    gd.addChoice("Sampling_method", samplingMethodNames, samplingMethodNames[samplingMethodIndex]);
    gd.addSlider("Sampling_factor", 0.2, 4, perimeterSamplingFactor);
    gd.addMessage("Estimation options:");
    String[] thresholdMethodNames = SettingsManager.getNames((Object[]) FRC.ThresholdMethod.values());
    gd.addChoice("Threshold_method", thresholdMethodNames, thresholdMethodNames[thresholdMethodIndex]);
    String[] precisionMethodNames = SettingsManager.getNames((Object[]) PrecisionMethod.values());
    gd.addChoice("Precision_method", precisionMethodNames, precisionMethodNames[precisionMethodIndex]);
    gd.addNumericField("Precision_Mean", mean, 2, 6, "nm");
    gd.addNumericField("Precision_Sigma", sigma, 2, 6, "nm");
    gd.addCheckbox("Sample_decay", sampleDecay);
    gd.addCheckbox("LOESS_smoothing", loessSmoothing);
    gd.addCheckbox("Fit_precision", fitPrecision);
    gd.addSlider("MinQ", 0, 0.4, minQ);
    gd.addSlider("MaxQ", 0.1, 0.5, maxQ);
    gd.showDialog();
    if (gd.wasCanceled())
        return false;
    inputOption = ResultsManager.getInputSource(gd);
    if (!titles.isEmpty())
        chooseRoi = gd.getNextBoolean();
    //useSignal = gd.getNextBoolean();
    imageScaleIndex = gd.getNextChoiceIndex();
    imageSizeIndex = gd.getNextChoiceIndex();
    blockSize = Math.max(1, (int) gd.getNextNumber());
    randomSplit = gd.getNextBoolean();
    maxPerBin = Math.abs((int) gd.getNextNumber());
    fourierMethodIndex = gd.getNextChoiceIndex();
    fourierMethod = FourierMethod.values()[fourierMethodIndex];
    samplingMethodIndex = gd.getNextChoiceIndex();
    samplingMethod = SamplingMethod.values()[samplingMethodIndex];
    perimeterSamplingFactor = gd.getNextNumber();
    thresholdMethodIndex = gd.getNextChoiceIndex();
    thresholdMethod = FRC.ThresholdMethod.values()[thresholdMethodIndex];
    precisionMethodIndex = gd.getNextChoiceIndex();
    precisionMethod = PrecisionMethod.values()[precisionMethodIndex];
    mean = Math.abs(gd.getNextNumber());
    sigma = Math.abs(gd.getNextNumber());
    sampleDecay = gd.getNextBoolean();
    loessSmoothing = gd.getNextBoolean();
    fitPrecision = gd.getNextBoolean();
    minQ = Maths.clip(0, 0.5, gd.getNextNumber());
    maxQ = Maths.clip(0, 0.5, gd.getNextNumber());
    // Check arguments
    try {
        Parameters.isAboveZero("Perimeter sampling factor", perimeterSamplingFactor);
        if (precisionMethod == PrecisionMethod.FIXED) {
            Parameters.isAboveZero("Precision Mean", mean);
            Parameters.isAboveZero("Precision Sigma", sigma);
        }
        Parameters.isAbove("MaxQ", maxQ, minQ);
    } catch (IllegalArgumentException e) {
        IJ.error(TITLE, e.getMessage());
        return false;
    }
    if (!titles.isEmpty() && chooseRoi) {
        if (titles.size() == 1) {
            roiImage = titles.get(0);
            Recorder.recordOption("Image", roiImage);
        } else {
            String[] items = titles.toArray(new String[titles.size()]);
            gd = new ExtendedGenericDialog(TITLE);
            gd.addMessage("Select the source image for the ROI");
            gd.addChoice("Image", items, roiImage);
            gd.showDialog();
            if (gd.wasCanceled())
                return false;
            roiImage = gd.getNextChoice();
        }
        ImagePlus imp = WindowManager.getImage(roiImage);
        roiBounds = imp.getRoi().getBounds();
        roiImageWidth = imp.getWidth();
        roiImageHeight = imp.getHeight();
    } else {
        roiBounds = null;
    }
    return true;
}