Search in sources :

Example 1 with BufferedTextWindow

use of gdsc.core.ij.BufferedTextWindow in project GDSC-SMLM by aherbert.

the class BenchmarkFilterAnalysis method filterAnalysis.

private int filterAnalysis(FilterSet filterSet, int setNumber, DirectFilter currentOptimum, double rangeReduction) {
    // Check if the filters are the same so allowing optimisation
    final boolean allSameType = filterSet.allSameType();
    this.ga_resultsList = resultsList;
    Chromosome<FilterScore> best = null;
    String algorithm = "";
    // All the search algorithms use search dimensions.
    // Create search dimensions if needed (these are used for testing if the optimum is at the limit).
    ss_filter = null;
    ss_lower = null;
    ss_upper = null;
    FixedDimension[] originalDimensions = null;
    boolean rangeInput = false;
    boolean[] disabled = null;
    double[][] seed = null;
    boolean nonInteractive = false;
    if (allSameType) {
        // There should always be 1 filter
        ss_filter = (DirectFilter) filterSet.getFilters().get(0);
        int n = ss_filter.getNumberOfParameters();
        // Option to configure a range
        rangeInput = filterSet.getName().contains("Range");
        double[] range = new double[n];
        if (rangeInput && filterSet.size() == 4) {
            originalDimensions = new FixedDimension[n];
            // This is used as min/lower/upper/max
            final Filter minF = ss_filter;
            final Filter lowerF = filterSet.getFilters().get(1);
            final Filter upperF = filterSet.getFilters().get(2);
            final Filter maxF = filterSet.getFilters().get(3);
            for (int i = 0; i < n; i++) {
                double min = minF.getParameterValue(i);
                double lower = lowerF.getParameterValue(i);
                double upper = upperF.getParameterValue(i);
                range[i] = upper - lower;
                double max = maxF.getParameterValue(i);
                double minIncrement = ss_filter.getParameterIncrement(i);
                try {
                    originalDimensions[i] = new FixedDimension(min, max, minIncrement, lower, upper);
                } catch (IllegalArgumentException e) {
                    Utils.log(TITLE + " : Unable to configure dimension [%d] %s: " + e.getMessage(), i, ss_filter.getParameterName(i));
                    originalDimensions = null;
                    rangeInput = false;
                    break;
                }
            }
        }
        if (rangeInput && (filterSet.size() == 3 || filterSet.size() == 2)) {
            originalDimensions = new FixedDimension[n];
            // This is used as lower/upper[/increment]
            final Filter lowerF = ss_filter;
            final Filter upperF = filterSet.getFilters().get(1);
            for (int i = 0; i < n; i++) {
                // Do not disable if the increment is not set. This is left to the user to decide.
                //					if (incF.getParameterValue(i) == incF.getDisabledParameterValue(i) ||
                //							Double.isInfinite(incF.getParameterValue(i)))
                //					{
                //						// Not enabled
                //						dimensions[i] = new SearchDimension(incF.getDisabledParameterValue(i));
                //						continue;
                //					}
                double lower = lowerF.getParameterValue(i);
                double upper = upperF.getParameterValue(i);
                range[i] = upper - lower;
                ParameterType type = ss_filter.getParameterType(i);
                double min = BenchmarkSpotFit.getMin(type);
                double max = BenchmarkSpotFit.getMax(type);
                double minIncrement = ss_filter.getParameterIncrement(i);
                try {
                    originalDimensions[i] = new FixedDimension(min, max, minIncrement, lower, upper);
                } catch (IllegalArgumentException e) {
                    Utils.log(TITLE + " : Unable to configure dimension [%d] %s: " + e.getMessage(), i, ss_filter.getParameterName(i));
                    originalDimensions = null;
                    rangeInput = false;
                    break;
                }
            }
        }
        // Get the dimensions from the filters
        if (originalDimensions == null) {
            originalDimensions = new FixedDimension[n];
            // Allow inputing a filter set (e.g. saved from previous optimisation)
            // Find the limits in the current scores
            final double[] lower = ss_filter.getParameters().clone();
            final double[] upper = lower.clone();
            // Allow the SearchSpace algorithms to be seeded with an initial population 
            // for the first evaluation of the optimum. This is done when the input filter 
            // set is not a range.
            seed = new double[filterSet.size()][];
            int c = 0;
            for (Filter f : filterSet.getFilters()) {
                final double[] point = f.getParameters();
                seed[c++] = point;
                for (int j = 0; j < lower.length; j++) {
                    if (lower[j] > point[j])
                        lower[j] = point[j];
                    if (upper[j] < point[j])
                        upper[j] = point[j];
                }
            }
            // Min/max must be set using values from BenchmarkSpotFit.
            for (int i = 0; i < n; i++) {
                if (lower[i] == upper[i]) {
                    // Not enabled
                    originalDimensions[i] = new FixedDimension(lower[i]);
                    continue;
                }
                ParameterType type = ss_filter.getParameterType(i);
                double min = BenchmarkSpotFit.getMin(type);
                double max = BenchmarkSpotFit.getMax(type);
                double minIncrement = ss_filter.getParameterIncrement(i);
                if (min > lower[i])
                    min = lower[i];
                if (max < upper[i])
                    max = upper[i];
                try {
                    originalDimensions[i] = new FixedDimension(min, max, minIncrement, lower[i], upper[i]);
                } catch (IllegalArgumentException e) {
                    Utils.log(TITLE + " : Unable to configure dimension [%d] %s: " + e.getMessage(), i, ss_filter.getParameterName(i));
                    originalDimensions = null;
                    break;
                }
            }
            if (originalDimensions == null) {
                // Failed to work out the dimensions. No optimisation will be possible.
                // Sort so that the filters are in a nice order for reporting
                filterSet.sort();
                // This will not be used when the dimensions are null
                seed = null;
            }
        }
        if (originalDimensions != null) {
            // Use the current optimum if we are doing a range optimisation
            if (currentOptimum != null && rangeInput && currentOptimum.getType().equals(ss_filter.getType()) && evolve != 0) {
                // Suppress dialogs and use the current settings
                nonInteractive = true;
                double[] p = currentOptimum.getParameters();
                // Range search uses SearchDimension and we must centre on the optimum after creation.							
                for (int i = 0; i < originalDimensions.length; i++) {
                    double centre = p[i];
                    double r = 0;
                    if (originalDimensions[i].isActive()) {
                        // Set the range around the centre.
                        // This uses the range for each param when we read the filters.
                        r = range[i];
                        // Optionally reduce the width of the dimensions. 
                        if (rangeReduction > 0 && rangeReduction < 1)
                            r *= rangeReduction;
                    }
                    double lower = centre - r * 0.5;
                    double upper = centre + r * 0.5;
                    originalDimensions[i] = originalDimensions[i].create(lower, upper);
                }
            }
            // Store the dimensions so we can do an 'at limit' check
            disabled = new boolean[originalDimensions.length];
            ss_lower = new double[originalDimensions.length];
            ss_upper = new double[originalDimensions.length];
            for (int i = 0; i < disabled.length; i++) {
                disabled[i] = !originalDimensions[i].isActive();
                ss_lower[i] = originalDimensions[i].lower;
                ss_upper[i] = originalDimensions[i].upper;
            }
        }
    } else {
        // Sort so that the filters are in a nice order for reporting
        filterSet.sort();
    }
    analysisStopWatch = StopWatch.createStarted();
    if (evolve == 1 && originalDimensions != null) {
        // Collect parameters for the genetic algorithm
        pauseFilterTimer();
        // Remember the step size settings
        double[] stepSize = stepSizeMap.get(setNumber);
        if (stepSize == null || stepSize.length != ss_filter.length()) {
            stepSize = ss_filter.mutationStepRange().clone();
            for (int j = 0; j < stepSize.length; j++) stepSize[j] *= delta;
            // See if the same number of parameters have been optimised in other algorithms
            boolean[] enabled = searchRangeMap.get(setNumber);
            if (enabled != null && enabled.length == stepSize.length) {
                for (int j = 0; j < stepSize.length; j++) if (!enabled[j])
                    stepSize[j] *= -1;
            }
        }
        GenericDialog gd = null;
        int[] indices = ss_filter.getChromosomeParameters();
        boolean runAlgorithm = nonInteractive;
        if (!nonInteractive) {
            // Ask the user for the mutation step parameters.
            gd = new GenericDialog(TITLE);
            String prefix = setNumber + "_";
            gd.addMessage("Configure the genetic algorithm for [" + setNumber + "] " + filterSet.getName());
            gd.addNumericField(prefix + "Population_size", populationSize, 0);
            gd.addNumericField(prefix + "Failure_limit", failureLimit, 0);
            gd.addNumericField(prefix + "Tolerance", tolerance, -1);
            gd.addNumericField(prefix + "Converged_count", convergedCount, 0);
            gd.addSlider(prefix + "Mutation_rate", 0.05, 1, mutationRate);
            gd.addSlider(prefix + "Crossover_rate", 0.05, 1, crossoverRate);
            gd.addSlider(prefix + "Mean_children", 0.05, 3, meanChildren);
            gd.addSlider(prefix + "Selection_fraction", 0.05, 0.5, selectionFraction);
            gd.addCheckbox(prefix + "Ramped_selection", rampedSelection);
            gd.addCheckbox(prefix + "Save_option", saveOption);
            gd.addMessage("Configure the step size for each parameter");
            for (int j = 0; j < indices.length; j++) {
                // Do not mutate parameters that were not expanded, i.e. the input did not vary them.
                final double step = (originalDimensions[indices[j]].isActive()) ? stepSize[j] : 0;
                gd.addNumericField(getDialogName(prefix, ss_filter, indices[j]), step, 2);
            }
            gd.showDialog();
            runAlgorithm = !gd.wasCanceled();
        }
        if (runAlgorithm) {
            // Used to create random sample
            FixedDimension[] dimensions = Arrays.copyOf(originalDimensions, originalDimensions.length);
            if (!nonInteractive) {
                populationSize = (int) Math.abs(gd.getNextNumber());
                if (populationSize < 10)
                    populationSize = 10;
                failureLimit = (int) Math.abs(gd.getNextNumber());
                tolerance = gd.getNextNumber();
                // Allow negatives
                convergedCount = (int) gd.getNextNumber();
                mutationRate = Math.abs(gd.getNextNumber());
                crossoverRate = Math.abs(gd.getNextNumber());
                meanChildren = Math.abs(gd.getNextNumber());
                selectionFraction = Math.abs(gd.getNextNumber());
                rampedSelection = gd.getNextBoolean();
                saveOption = gd.getNextBoolean();
                for (int j = 0; j < indices.length; j++) {
                    stepSize[j] = gd.getNextNumber();
                }
                // Store for repeat analysis
                stepSizeMap.put(setNumber, stepSize);
            }
            for (int j = 0; j < indices.length; j++) {
                // A zero step size will keep the parameter but prevent range mutation.
                if (stepSize[j] < 0) {
                    dimensions[indices[j]] = new FixedDimension(ss_filter.getDisabledParameterValue(indices[j]));
                    disabled[indices[j]] = true;
                }
            }
            //				// Reset negatives to zero
            //				stepSize = stepSize.clone();
            //				for (int j = 0; j < stepSize.length; j++)
            //					if (stepSize[j] < 0)
            //						stepSize[j] = 0;
            // Create the genetic algorithm
            RandomDataGenerator random = new RandomDataGenerator(new Well44497b());
            SimpleMutator<FilterScore> mutator = new SimpleMutator<FilterScore>(random, mutationRate);
            // Override the settings with the step length, a min of zero and the configured upper
            double[] upper = ss_filter.upperLimit();
            mutator.overrideChromosomeSettings(stepSize, new double[stepSize.length], upper);
            Recombiner<FilterScore> recombiner = new SimpleRecombiner<FilterScore>(random, crossoverRate, meanChildren);
            SelectionStrategy<FilterScore> selectionStrategy;
            // If the initial population is huge ensure that the first selection culls to the correct size
            final int selectionMax = (int) (selectionFraction * populationSize);
            if (rampedSelection)
                selectionStrategy = new RampedSelectionStrategy<FilterScore>(random, selectionFraction, selectionMax);
            else
                selectionStrategy = new SimpleSelectionStrategy<FilterScore>(random, selectionFraction, selectionMax);
            ToleranceChecker<FilterScore> ga_checker = new InterruptChecker(tolerance, tolerance * 1e-3, convergedCount);
            // Create new random filters if the population is initially below the population size
            List<Filter> filters = filterSet.getFilters();
            if (filterSet.size() < populationSize) {
                filters = new ArrayList<Filter>(populationSize);
                // Add the existing filters if they are not a range input file
                if (!rangeInput)
                    filters.addAll(filterSet.getFilters());
                // Add current optimum to seed
                if (nonInteractive)
                    filters.add(currentOptimum);
                // The GA does not use the min interval grid so sample without rounding
                double[][] sample = SearchSpace.sampleWithoutRounding(dimensions, populationSize - filters.size(), null);
                filters.addAll(searchSpaceToFilters(sample));
            }
            ga_population = new Population<FilterScore>(filters);
            ga_population.setPopulationSize(populationSize);
            ga_population.setFailureLimit(failureLimit);
            selectionStrategy.setTracker(this);
            // Evolve
            algorithm = EVOLVE[evolve];
            ga_statusPrefix = algorithm + " [" + setNumber + "] " + filterSet.getName() + " ... ";
            ga_iteration = 0;
            ga_population.setTracker(this);
            createGAWindow();
            resumeFilterTimer();
            best = ga_population.evolve(mutator, recombiner, this, selectionStrategy, ga_checker);
            if (best != null) {
                // In case optimisation was stopped
                IJ.resetEscape();
                // The GA may produce coordinates off the min interval grid
                best = enumerateMinInterval(best, stepSize, indices);
                // Now update the filter set for final assessment
                filterSet = new FilterSet(filterSet.getName(), populationToFilters(ga_population.getIndividuals()));
                // Option to save the filters
                if (saveOption)
                    saveFilterSet(filterSet, setNumber, !nonInteractive);
            }
        } else
            resumeFilterTimer();
    }
    if ((evolve == 2 || evolve == 4) && originalDimensions != null) {
        // Collect parameters for the range search algorithm
        pauseFilterTimer();
        boolean isStepSearch = evolve == 4;
        // The step search should use a multi-dimension refinement and no range reduction
        SearchSpace.RefinementMode myRefinementMode = SearchSpace.RefinementMode.MULTI_DIMENSION;
        // Remember the enabled settings
        boolean[] enabled = searchRangeMap.get(setNumber);
        int n = ss_filter.getNumberOfParameters();
        if (enabled == null || enabled.length != n) {
            enabled = new boolean[n];
            Arrays.fill(enabled, true);
            // See if the same number of parameters have been optimised in other algorithms
            double[] stepSize = stepSizeMap.get(setNumber);
            if (stepSize != null && enabled.length == stepSize.length) {
                for (int j = 0; j < stepSize.length; j++) if (stepSize[j] < 0)
                    enabled[j] = false;
            }
        }
        GenericDialog gd = null;
        boolean runAlgorithm = nonInteractive;
        if (!nonInteractive) {
            // Ask the user for the search parameters.
            gd = new GenericDialog(TITLE);
            String prefix = setNumber + "_";
            gd.addMessage("Configure the " + EVOLVE[evolve] + " algorithm for [" + setNumber + "] " + filterSet.getName());
            gd.addSlider(prefix + "Width", 1, 5, rangeSearchWidth);
            if (!isStepSearch) {
                gd.addCheckbox(prefix + "Save_option", saveOption);
                gd.addNumericField(prefix + "Max_iterations", maxIterations, 0);
                String[] modes = SettingsManager.getNames((Object[]) SearchSpace.RefinementMode.values());
                gd.addSlider(prefix + "Reduce", 0.01, 0.99, rangeSearchReduce);
                gd.addChoice("Refinement", modes, modes[refinementMode]);
            }
            gd.addNumericField(prefix + "Seed_size", seedSize, 0);
            // Add choice of fields to optimise
            for (int i = 0; i < n; i++) gd.addCheckbox(getDialogName(prefix, ss_filter, i), enabled[i]);
            gd.showDialog();
            runAlgorithm = !gd.wasCanceled();
        }
        if (runAlgorithm) {
            SearchDimension[] dimensions = new SearchDimension[n];
            if (!nonInteractive) {
                rangeSearchWidth = (int) gd.getNextNumber();
                if (!isStepSearch) {
                    saveOption = gd.getNextBoolean();
                    maxIterations = (int) gd.getNextNumber();
                    refinementMode = gd.getNextChoiceIndex();
                    rangeSearchReduce = gd.getNextNumber();
                }
                seedSize = (int) gd.getNextNumber();
                for (int i = 0; i < n; i++) enabled[i] = gd.getNextBoolean();
                searchRangeMap.put(setNumber, enabled);
            }
            if (!isStepSearch)
                myRefinementMode = SearchSpace.RefinementMode.values()[refinementMode];
            for (int i = 0; i < n; i++) {
                if (enabled[i]) {
                    try {
                        dimensions[i] = originalDimensions[i].create(rangeSearchWidth);
                        dimensions[i].setPad(true);
                        // Prevent range reduction so that the step search just does a single refinement step
                        dimensions[i].setReduceFactor((isStepSearch) ? 1 : rangeSearchReduce);
                        // Centre on current optimum
                        if (nonInteractive)
                            dimensions[i].setCentre(currentOptimum.getParameterValue(i));
                    } catch (IllegalArgumentException e) {
                        IJ.error(TITLE, String.format("Unable to configure dimension [%d] %s: " + e.getMessage(), i, ss_filter.getParameterName(i)));
                        return -1;
                    }
                } else {
                    dimensions[i] = new SearchDimension(ss_filter.getDisabledParameterValue(i));
                }
            }
            for (int i = 0; i < disabled.length; i++) disabled[i] = !dimensions[i].active;
            // Check the number of combinations is OK
            long combinations = SearchSpace.countCombinations(dimensions);
            if (!nonInteractive && combinations > 10000) {
                gd = new GenericDialog(TITLE);
                gd.addMessage(String.format("%d combinations for the configured dimensions.\n \nClick 'Yes' to optimise.", combinations));
                gd.enableYesNoCancel();
                gd.hideCancelButton();
                gd.showDialog();
                if (!gd.wasOKed()) {
                    combinations = 0;
                }
            }
            if (combinations == 0) {
                resumeFilterTimer();
            } else {
                algorithm = EVOLVE[evolve] + " " + rangeSearchWidth;
                ga_statusPrefix = algorithm + " [" + setNumber + "] " + filterSet.getName() + " ... ";
                ga_iteration = 0;
                es_optimum = null;
                SearchSpace ss = new SearchSpace();
                ss.setTracker(this);
                if (seedSize > 0) {
                    double[][] sample;
                    // Add current optimum to seed
                    if (nonInteractive) {
                        sample = new double[1][];
                        sample[0] = currentOptimum.getParameters();
                        seed = merge(seed, sample);
                    }
                    int size = (seed == null) ? 0 : seed.length;
                    // Sample without rounding as the seed will be rounded
                    sample = SearchSpace.sampleWithoutRounding(dimensions, seedSize - size, null);
                    seed = merge(seed, sample);
                }
                // Note: If we have an optimum and we are not seeding this should not matter as the dimensions 
                // have been centred on the current optimum					
                ss.seed(seed);
                ConvergenceChecker<FilterScore> checker = new InterruptConvergenceChecker(0, 0, maxIterations);
                createGAWindow();
                resumeFilterTimer();
                SearchResult<FilterScore> optimum = ss.search(dimensions, this, checker, myRefinementMode);
                if (optimum != null) {
                    // In case optimisation was stopped
                    IJ.resetEscape();
                    best = ((SimpleFilterScore) optimum.score).r.filter;
                    if (seedSize > 0) {
                    // Not required as the search now respects the min interval
                    // The optimum may be off grid if it was from the seed
                    //best = enumerateMinInterval(best, enabled);
                    }
                    // Now update the filter set for final assessment
                    filterSet = new FilterSet(filterSet.getName(), searchSpaceToFilters((DirectFilter) best, ss.getSearchSpace()));
                    // Option to save the filters
                    if (saveOption)
                        saveFilterSet(filterSet, setNumber, !nonInteractive);
                }
            }
        } else
            resumeFilterTimer();
    }
    if (evolve == 3 && originalDimensions != null) {
        // Collect parameters for the enrichment search algorithm
        pauseFilterTimer();
        boolean[] enabled = searchRangeMap.get(setNumber);
        int n = ss_filter.getNumberOfParameters();
        if (enabled == null || enabled.length != n) {
            enabled = new boolean[n];
            Arrays.fill(enabled, true);
            // See if the same number of parameters have been optimised in other algorithms
            double[] stepSize = stepSizeMap.get(setNumber);
            if (stepSize != null && enabled.length == stepSize.length) {
                for (int j = 0; j < stepSize.length; j++) if (stepSize[j] < 0)
                    enabled[j] = false;
            }
        }
        GenericDialog gd = null;
        boolean runAlgorithm = nonInteractive;
        if (!nonInteractive) {
            // Ask the user for the search parameters.
            gd = new GenericDialog(TITLE);
            String prefix = setNumber + "_";
            gd.addMessage("Configure the enrichment search algorithm for [" + setNumber + "] " + filterSet.getName());
            gd.addCheckbox(prefix + "Save_option", saveOption);
            gd.addNumericField(prefix + "Max_iterations", maxIterations, 0);
            gd.addNumericField(prefix + "Converged_count", convergedCount, 0);
            gd.addNumericField(prefix + "Samples", enrichmentSamples, 0);
            gd.addSlider(prefix + "Fraction", 0.01, 0.99, enrichmentFraction);
            gd.addSlider(prefix + "Padding", 0, 0.99, enrichmentPadding);
            // Add choice of fields to optimise
            for (int i = 0; i < n; i++) gd.addCheckbox(getDialogName(prefix, ss_filter, i), enabled[i]);
            gd.showDialog();
            runAlgorithm = !gd.wasCanceled();
        }
        if (runAlgorithm) {
            FixedDimension[] dimensions = Arrays.copyOf(originalDimensions, originalDimensions.length);
            if (!nonInteractive) {
                saveOption = gd.getNextBoolean();
                maxIterations = (int) gd.getNextNumber();
                convergedCount = (int) gd.getNextNumber();
                enrichmentSamples = (int) gd.getNextNumber();
                enrichmentFraction = gd.getNextNumber();
                enrichmentPadding = gd.getNextNumber();
                for (int i = 0; i < n; i++) enabled[i] = gd.getNextBoolean();
                searchRangeMap.put(setNumber, enabled);
            }
            for (int i = 0; i < n; i++) {
                if (!enabled[i])
                    dimensions[i] = new FixedDimension(ss_filter.getDisabledParameterValue(i));
            }
            for (int i = 0; i < disabled.length; i++) disabled[i] = !dimensions[i].active;
            algorithm = EVOLVE[evolve];
            ga_statusPrefix = algorithm + " [" + setNumber + "] " + filterSet.getName() + " ... ";
            ga_iteration = 0;
            es_optimum = null;
            SearchSpace ss = new SearchSpace();
            ss.setTracker(this);
            // Add current optimum to seed
            if (nonInteractive) {
                double[][] sample = new double[1][];
                sample[0] = currentOptimum.getParameters();
                seed = merge(seed, sample);
            }
            ss.seed(seed);
            ConvergenceChecker<FilterScore> checker = new InterruptConvergenceChecker(0, 0, maxIterations, convergedCount);
            createGAWindow();
            resumeFilterTimer();
            SearchResult<FilterScore> optimum = ss.enrichmentSearch(dimensions, this, checker, enrichmentSamples, enrichmentFraction, enrichmentPadding);
            if (optimum != null) {
                // In case optimisation was stopped
                IJ.resetEscape();
                best = ((SimpleFilterScore) optimum.score).r.filter;
                // Not required as the search now respects the min interval
                // Enumerate on the min interval to produce the final filter
                //best = enumerateMinInterval(best, enabled);
                // Now update the filter set for final assessment
                filterSet = new FilterSet(filterSet.getName(), searchSpaceToFilters((DirectFilter) best, ss.getSearchSpace()));
                // Option to save the filters
                if (saveOption)
                    saveFilterSet(filterSet, setNumber, !nonInteractive);
            }
        } else
            resumeFilterTimer();
    }
    IJ.showStatus("Analysing [" + setNumber + "] " + filterSet.getName() + " ...");
    // Do not support plotting if we used optimisation
    double[] xValues = (best != null || isHeadless || (plotTopN == 0)) ? null : new double[filterSet.size()];
    double[] yValues = (xValues == null) ? null : new double[xValues.length];
    SimpleFilterScore max = null;
    // It can just assess the top 1 required for the summary.
    if (best != null) {
        // Only assess the top 1 filter for the summary
        List<Filter> list = new ArrayList<Filter>();
        list.add((DirectFilter) best);
        filterSet = new FilterSet(filterSet.getName(), list);
    }
    // Score the filters and report the results if configured.
    FilterScoreResult[] scoreResults = scoreFilters(setUncomputedStrength(filterSet), showResultsTable);
    if (scoreResults == null)
        return -1;
    analysisStopWatch.stop();
    for (int index = 0; index < scoreResults.length; index++) {
        final FilterScoreResult scoreResult = scoreResults[index];
        if (xValues != null) {
            xValues[index] = scoreResult.filter.getNumericalValue();
            yValues[index] = scoreResult.score;
        }
        final SimpleFilterScore result = new SimpleFilterScore(scoreResult, allSameType, scoreResult.criteria >= minCriteria);
        if (result.compareTo(max) < 0) {
            max = result;
        }
    }
    if (showResultsTable) {
        BufferedTextWindow tw = null;
        if (resultsWindow != null) {
            tw = new BufferedTextWindow(resultsWindow);
            tw.setIncrement(Integer.MAX_VALUE);
        }
        for (int index = 0; index < scoreResults.length; index++) addToResultsWindow(tw, scoreResults[index].text);
        if (resultsWindow != null)
            resultsWindow.getTextPanel().updateDisplay();
    }
    // Check the top filter against the limits of the original dimensions
    char[] atLimit = null;
    if (allSameType && originalDimensions != null) {
        DirectFilter filter = max.r.filter;
        int[] indices = filter.getChromosomeParameters();
        atLimit = new char[indices.length];
        StringBuilder sb = new StringBuilder(200);
        for (int j = 0; j < indices.length; j++) {
            atLimit[j] = ComplexFilterScore.WITHIN;
            final int p = indices[j];
            if (disabled[p])
                continue;
            final double value = filter.getParameterValue(p);
            double lowerLimit = originalDimensions[p].getLower();
            double upperLimit = originalDimensions[p].getUpper();
            int c1 = Double.compare(value, lowerLimit);
            if (c1 <= 0) {
                atLimit[j] = ComplexFilterScore.FLOOR;
                sb.append(" : ").append(filter.getParameterName(p)).append(' ').append(atLimit[j]).append('[').append(Utils.rounded(value));
                if (c1 == -1) {
                    atLimit[j] = ComplexFilterScore.BELOW;
                    sb.append("<").append(Utils.rounded(lowerLimit));
                }
                sb.append("]");
            } else {
                int c2 = Double.compare(value, upperLimit);
                if (c2 >= 0) {
                    atLimit[j] = ComplexFilterScore.CEIL;
                    sb.append(" : ").append(filter.getParameterName(p)).append(' ').append(atLimit[j]).append('[').append(Utils.rounded(value));
                    if (c2 == 1) {
                        atLimit[j] = ComplexFilterScore.ABOVE;
                        sb.append(">").append(Utils.rounded(upperLimit));
                    }
                    sb.append("]");
                }
            }
        }
        if (sb.length() > 0) {
            if (max.criteriaPassed) {
                Utils.log("Warning: Top filter (%s @ %s|%s) [%s] at the limit of the expanded range%s", filter.getName(), Utils.rounded((invertScore) ? -max.score : max.score), Utils.rounded((invertCriteria) ? -minCriteria : minCriteria), limitFailCount + limitRange, sb.toString());
            } else {
                Utils.log("Warning: Top filter (%s @ -|%s) [%s] at the limit of the expanded range%s", filter.getName(), Utils.rounded((invertCriteria) ? -max.criteria : max.criteria), limitFailCount + limitRange, sb.toString());
            }
        }
    }
    // Note that max should never be null since this method is not run with an empty filter set
    // We may have no filters that pass the criteria
    String type = max.r.filter.getType();
    if (!max.criteriaPassed) {
        Utils.log("Warning: Filter does not pass the criteria: %s : Best = %s using %s", type, Utils.rounded((invertCriteria) ? -max.criteria : max.criteria), max.r.filter.getName());
        return 0;
    }
    // This could be an option?
    boolean allowDuplicates = true;
    // XXX - Commented out the requirement to be the same type to store for later analysis. 
    // This may break the code, however I think that all filter sets should be able to have a best filter
    // irrespective of whether they were the same type or not.
    //if (allSameType)
    //{
    ComplexFilterScore newFilterScore = new ComplexFilterScore(max.r, atLimit, algorithm, analysisStopWatch.getTime(), "", 0);
    addBestFilter(type, allowDuplicates, newFilterScore);
    // Add spacer at end of each result set
    if (isHeadless) {
        if (showResultsTable && filterSet.size() > 1)
            IJ.log("");
    } else {
        if (showResultsTable && filterSet.size() > 1)
            resultsWindow.append("");
        if (plotTopN > 0 && xValues != null) {
            // Check the xValues are unique. Since the filters have been sorted by their
            // numeric value we only need to compare adjacent entries.
            boolean unique = true;
            for (int ii = 0; ii < xValues.length - 1; ii++) {
                if (xValues[ii] == xValues[ii + 1]) {
                    unique = false;
                    break;
                }
            }
            String xAxisName = filterSet.getValueName();
            if (unique) {
                // Check the values all refer to the same property
                for (Filter filter : filterSet.getFilters()) {
                    if (!xAxisName.equals(filter.getNumericalValueName())) {
                        unique = false;
                        break;
                    }
                }
            }
            if (!unique) {
                // If not unique then renumber them and use an arbitrary label
                xAxisName = "Filter";
                for (int ii = 0; ii < xValues.length; ii++) xValues[ii] = ii + 1;
            }
            String title = filterSet.getName();
            // Check if a previous filter set had the same name, update if necessary
            NamedPlot p = getNamedPlot(title);
            if (p == null)
                plots.add(new NamedPlot(title, xAxisName, xValues, yValues));
            else
                p.updateValues(xAxisName, xValues, yValues);
            if (plots.size() > plotTopN) {
                Collections.sort(plots);
                p = plots.remove(plots.size() - 1);
            }
        }
    }
    return 0;
}
Also used : SimpleRecombiner(gdsc.smlm.ga.SimpleRecombiner) SearchSpace(gdsc.smlm.search.SearchSpace) ArrayList(java.util.ArrayList) RandomDataGenerator(org.apache.commons.math3.random.RandomDataGenerator) BufferedTextWindow(gdsc.core.ij.BufferedTextWindow) SearchDimension(gdsc.smlm.search.SearchDimension) FixedDimension(gdsc.smlm.search.FixedDimension) FilterSet(gdsc.smlm.results.filter.FilterSet) RampedSelectionStrategy(gdsc.smlm.ga.RampedSelectionStrategy) GenericDialog(ij.gui.GenericDialog) NonBlockingGenericDialog(ij.gui.NonBlockingGenericDialog) ParameterType(gdsc.smlm.results.filter.ParameterType) IDirectFilter(gdsc.smlm.results.filter.IDirectFilter) DirectFilter(gdsc.smlm.results.filter.DirectFilter) SimpleSelectionStrategy(gdsc.smlm.ga.SimpleSelectionStrategy) IDirectFilter(gdsc.smlm.results.filter.IDirectFilter) DirectFilter(gdsc.smlm.results.filter.DirectFilter) Filter(gdsc.smlm.results.filter.Filter) MultiPathFilter(gdsc.smlm.results.filter.MultiPathFilter) MaximaSpotFilter(gdsc.smlm.filters.MaximaSpotFilter) SimpleMutator(gdsc.smlm.ga.SimpleMutator) FilterScore(gdsc.smlm.results.filter.FilterScore) Well44497b(org.apache.commons.math3.random.Well44497b)

Example 2 with BufferedTextWindow

use of gdsc.core.ij.BufferedTextWindow in project GDSC-SMLM by aherbert.

the class BackgroundEstimator method showResults.

private void showResults() {
    Collections.sort(results, new Comparator<double[]>() {

        public int compare(double[] o1, double[] o2) {
            // Sort on slice number
            return (o1[0] < o2[0]) ? -1 : 1;
        }
    });
    BufferedTextWindow tw = new BufferedTextWindow(new TextWindow(imp.getTitle() + " Background", createHeader(), "", 800, 400));
    for (double[] result : results) tw.append(createResult(result));
    tw.flush();
}
Also used : BufferedTextWindow(gdsc.core.ij.BufferedTextWindow) TextWindow(ij.text.TextWindow) BufferedTextWindow(gdsc.core.ij.BufferedTextWindow)

Example 3 with BufferedTextWindow

use of gdsc.core.ij.BufferedTextWindow in project GDSC-SMLM by aherbert.

the class BenchmarkSpotFilter method getTable.

private BufferedTextWindow getTable(boolean batchSummary) {
    if (batchSummary) {
        if (batchSummaryTable == null || !batchSummaryTable.isVisible()) {
            TextWindow table = new TextWindow(TITLE + " Batch", createHeader(), "", 1000, 300);
            table.setVisible(true);
            batchSummaryTable = new BufferedTextWindow(table);
        }
        return batchSummaryTable;
    } else {
        if (summaryTable == null || !summaryTable.isVisible()) {
            TextWindow table = new TextWindow(TITLE, createHeader(), "", 1000, 300);
            table.setVisible(true);
            summaryTable = new BufferedTextWindow(table);
        }
        return summaryTable;
    }
}
Also used : TextWindow(ij.text.TextWindow) BufferedTextWindow(gdsc.core.ij.BufferedTextWindow) BufferedTextWindow(gdsc.core.ij.BufferedTextWindow)

Example 4 with BufferedTextWindow

use of gdsc.core.ij.BufferedTextWindow in project GDSC-SMLM by aherbert.

the class BenchmarkSpotFilter method summariseResults.

private BenchmarkFilterResult summariseResults(TIntObjectHashMap<FilterResult> filterResults, FitEngineConfiguration config, MaximaSpotFilter spotFilter, boolean relativeDistances, boolean batchSummary) {
    BenchmarkFilterResult filterResult = new BenchmarkFilterResult(filterResults, config, spotFilter);
    // Note: 
    // Although we can compute the TP/FP score as each additional spot is added
    // using the RankedScoreCalculator this is not applicable to the PeakFit method.
    // The method relies on all spot candidates being present in order to make a
    // decision to fit the candidate as a multiple. So scoring the filter candidates using
    // for example the top 10 may get a better score than if all candidates were scored
    // and the scores accumulated for the top 10, it is not how the algorithm will use the 
    // candidate set. I.e. It does not use the top 10, then top 20 to refine the fit, etc. 
    // (the method is not iterative) .
    // We require an assessment of how a subset of the scored candidates
    // in ranked order contributes to the overall score, i.e. are the candidates ranked
    // in the correct order, those most contributing to the match to the underlying data 
    // should be higher up and those least contributing will be at the end.
    // TODO We could add some smart filtering of candidates before ranking. This would
    // allow assessment of the candidate set handed to PeakFit. E.g. Threshold the image
    // and only use candidates that are in the foreground region.
    double[][] cumul = histogramFailures(filterResult);
    // Create the overall match score
    final double[] total = new double[3];
    final ArrayList<ScoredSpot> allSpots = new ArrayList<BenchmarkSpotFilter.ScoredSpot>();
    filterResults.forEachValue(new TObjectProcedure<FilterResult>() {

        public boolean execute(FilterResult result) {
            total[0] += result.result.getTP();
            total[1] += result.result.getFP();
            total[2] += result.result.getFN();
            allSpots.addAll(Arrays.asList(result.spots));
            return true;
        }
    });
    double tp = total[0], fp = total[1], fn = total[2];
    FractionClassificationResult allResult = new FractionClassificationResult(tp, fp, 0, fn);
    // The number of actual results
    final double n = (tp + fn);
    StringBuilder sb = new StringBuilder();
    double signal = (simulationParameters.minSignal + simulationParameters.maxSignal) * 0.5;
    // Create the benchmark settings and the fitting settings
    sb.append(imp.getStackSize()).append("\t");
    final int w = lastAnalysisBorder.width;
    final int h = lastAnalysisBorder.height;
    sb.append(w).append("\t");
    sb.append(h).append("\t");
    sb.append(Utils.rounded(n)).append("\t");
    double density = (n / imp.getStackSize()) / (w * h) / (simulationParameters.a * simulationParameters.a / 1e6);
    sb.append(Utils.rounded(density)).append("\t");
    sb.append(Utils.rounded(signal)).append("\t");
    sb.append(Utils.rounded(simulationParameters.s)).append("\t");
    sb.append(Utils.rounded(simulationParameters.a)).append("\t");
    sb.append(Utils.rounded(simulationParameters.depth)).append("\t");
    sb.append(simulationParameters.fixedDepth).append("\t");
    sb.append(Utils.rounded(simulationParameters.gain)).append("\t");
    sb.append(Utils.rounded(simulationParameters.readNoise)).append("\t");
    sb.append(Utils.rounded(simulationParameters.b)).append("\t");
    sb.append(Utils.rounded(simulationParameters.b2)).append("\t");
    // Compute the noise
    double noise = simulationParameters.b2;
    if (simulationParameters.emCCD) {
        // The b2 parameter was computed without application of the EM-CCD noise factor of 2.
        //final double b2 = backgroundVariance + readVariance
        //                = simulationParameters.b + readVariance
        // This should be applied only to the background variance.
        final double readVariance = noise - simulationParameters.b;
        noise = simulationParameters.b * 2 + readVariance;
    }
    sb.append(Utils.rounded(signal / Math.sqrt(noise))).append("\t");
    sb.append(Utils.rounded(simulationParameters.s / simulationParameters.a)).append("\t");
    sb.append(config.getDataFilterType()).append("\t");
    //sb.append(spotFilter.getName()).append("\t");
    sb.append(spotFilter.getSearch()).append("\t");
    sb.append(spotFilter.getBorder()).append("\t");
    sb.append(Utils.rounded(spotFilter.getSpread())).append("\t");
    sb.append(config.getDataFilter(0)).append("\t");
    final double param = config.getSmooth(0);
    final double hwhmMin = config.getHWHMMin();
    if (relativeDistances) {
        sb.append(Utils.rounded(param * hwhmMin)).append("\t");
        sb.append(Utils.rounded(param)).append("\t");
    } else {
        sb.append(Utils.rounded(param)).append("\t");
        sb.append(Utils.rounded(param / hwhmMin)).append("\t");
    }
    sb.append(spotFilter.getDescription()).append("\t");
    sb.append(lastAnalysisBorder.x).append("\t");
    sb.append(MATCHING_METHOD[matchingMethod]).append("\t");
    sb.append(Utils.rounded(lowerMatchDistance)).append("\t");
    sb.append(Utils.rounded(matchDistance)).append("\t");
    sb.append(Utils.rounded(lowerSignalFactor)).append("\t");
    sb.append(Utils.rounded(upperSignalFactor));
    resultPrefix = sb.toString();
    // Add the results
    sb.append("\t");
    // Rank the scored spots by intensity
    Collections.sort(allSpots);
    // Produce Recall, Precision, Jaccard for each cut of the spot candidates
    double[] r = new double[allSpots.size() + 1];
    double[] p = new double[r.length];
    double[] j = new double[r.length];
    double[] c = new double[r.length];
    double[] truePositives = new double[r.length];
    double[] falsePositives = new double[r.length];
    double[] intensity = new double[r.length];
    // Note: fn = n - tp
    tp = fp = 0;
    int i = 1;
    p[0] = 1;
    FastCorrelator corr = new FastCorrelator();
    double lastC = 0;
    double[] i1 = new double[r.length];
    double[] i2 = new double[r.length];
    int ci = 0;
    SimpleRegression regression = new SimpleRegression(false);
    for (ScoredSpot s : allSpots) {
        if (s.match) {
            // Score partial matches as part true-positive and part false-positive.
            // TP can be above 1 if we are allowing multiple matches.
            tp += s.getScore();
            fp += s.antiScore();
            // Just use a rounded intensity for now
            final double spotIntensity = s.getIntensity();
            final long v1 = (long) Math.round(spotIntensity);
            final long v2 = (long) Math.round(s.intensity);
            regression.addData(spotIntensity, s.intensity);
            i1[ci] = spotIntensity;
            i2[ci] = s.intensity;
            ci++;
            corr.add(v1, v2);
            lastC = corr.getCorrelation();
        } else
            fp++;
        r[i] = (double) tp / n;
        p[i] = (double) tp / (tp + fp);
        // (tp+fp+fn) == (fp+n) since tp+fn=n;
        j[i] = (double) tp / (fp + n);
        c[i] = lastC;
        truePositives[i] = tp;
        falsePositives[i] = fp;
        intensity[i] = s.getIntensity();
        i++;
    }
    i1 = Arrays.copyOf(i1, ci);
    i2 = Arrays.copyOf(i2, ci);
    final double slope = regression.getSlope();
    sb.append(Utils.rounded(slope)).append("\t");
    addResult(sb, allResult, c[c.length - 1]);
    // Output the match results when the recall achieves the fraction of the maximum.
    double target = r[r.length - 1];
    if (recallFraction < 100)
        target *= recallFraction / 100.0;
    int fractionIndex = 0;
    while (fractionIndex < r.length && r[fractionIndex] < target) {
        fractionIndex++;
    }
    if (fractionIndex == r.length)
        fractionIndex--;
    addResult(sb, new FractionClassificationResult(truePositives[fractionIndex], falsePositives[fractionIndex], 0, n - truePositives[fractionIndex]), c[fractionIndex]);
    // Output the match results at the maximum jaccard score
    int maxIndex = 0;
    for (int ii = 1; ii < r.length; ii++) {
        if (j[maxIndex] < j[ii])
            maxIndex = ii;
    }
    addResult(sb, new FractionClassificationResult(truePositives[maxIndex], falsePositives[maxIndex], 0, n - truePositives[maxIndex]), c[maxIndex]);
    sb.append(Utils.rounded(time / 1e6));
    // Calculate AUC (Average precision == Area Under Precision-Recall curve)
    final double auc = AUCCalculator.auc(p, r);
    // Compute the AUC using the adjusted precision curve
    // which uses the maximum precision for recall >= r
    final double[] maxp = new double[p.length];
    double max = 0;
    for (int k = maxp.length; k-- > 0; ) {
        if (max < p[k])
            max = p[k];
        maxp[k] = max;
    }
    final double auc2 = AUCCalculator.auc(maxp, r);
    sb.append("\t").append(Utils.rounded(auc));
    sb.append("\t").append(Utils.rounded(auc2));
    // Output the number of fit failures that must be processed to capture fractions of the true positives
    if (cumul[0].length != 0) {
        sb.append("\t").append(Utils.rounded(getFailures(cumul, 0.80)));
        sb.append("\t").append(Utils.rounded(getFailures(cumul, 0.90)));
        sb.append("\t").append(Utils.rounded(getFailures(cumul, 0.95)));
        sb.append("\t").append(Utils.rounded(getFailures(cumul, 0.99)));
        sb.append("\t").append(Utils.rounded(cumul[0][cumul[0].length - 1]));
    } else
        sb.append("\t\t\t\t\t");
    BufferedTextWindow resultsTable = getTable(batchSummary);
    resultsTable.append(sb.toString());
    // Store results
    filterResult.auc = auc;
    filterResult.auc2 = auc2;
    filterResult.r = r;
    filterResult.p = p;
    filterResult.j = j;
    filterResult.c = c;
    filterResult.maxIndex = maxIndex;
    filterResult.fractionIndex = fractionIndex;
    filterResult.cumul = cumul;
    filterResult.slope = slope;
    filterResult.i1 = i1;
    filterResult.i2 = i2;
    filterResult.intensity = intensity;
    filterResult.relativeDistances = relativeDistances;
    filterResult.time = time;
    return filterResult;
}
Also used : BufferedTextWindow(gdsc.core.ij.BufferedTextWindow) FastCorrelator(gdsc.core.utils.FastCorrelator) ArrayList(java.util.ArrayList) PeakResultPoint(gdsc.smlm.ij.plugins.ResultsMatchCalculator.PeakResultPoint) BasePoint(gdsc.core.match.BasePoint) SimpleRegression(org.apache.commons.math3.stat.regression.SimpleRegression) FractionClassificationResult(gdsc.core.match.FractionClassificationResult)

Example 5 with BufferedTextWindow

use of gdsc.core.ij.BufferedTextWindow in project GDSC-SMLM by aherbert.

the class BenchmarkFilterAnalysis method parameterAnalysis.

/**
	 * Run the optimum filter on a set of labelled peak results using various parameter settings outputting performance
	 * statistics on the success of the filter to an ImageJ table.
	 * <p>
	 * If a new optimum is found the class level static parameters are updated.
	 *
	 * @param nonInteractive
	 *            True if non interactive
	 * @param currentOptimum
	 *            the optimum
	 * @param rangeReduction
	 *            the range reduction
	 * @return the best filter
	 */
private ComplexFilterScore parameterAnalysis(boolean nonInteractive, ComplexFilterScore currentOptimum, double rangeReduction) {
    this.ga_resultsList = resultsList;
    String algorithm = "";
    // All the search algorithms use search dimensions.
    ss_filter = currentOptimum.r.filter;
    FixedDimension[] originalDimensions = new FixedDimension[3];
    double[] point = createParameters();
    String[] names = { "Fail count", "Residuals threshold", "Duplicate distance" };
    {
        // Local scope for i
        int i = 0;
        try {
            originalDimensions[i++] = new FixedDimension(minFailCount, maxFailCount, 1);
            // TODO - let the min intervals be configured, maybe via extra options
            if (BenchmarkSpotFit.computeDoublets)
                originalDimensions[i++] = new FixedDimension(minResidualsThreshold, maxResidualsThreshold, 0.05);
            else
                originalDimensions[i++] = new FixedDimension(1, 1, 0.05);
            originalDimensions[i++] = new FixedDimension(minDuplicateDistance, maxDuplicateDistance, 0.5);
        } catch (IllegalArgumentException e) {
            Utils.log(TITLE + " : Unable to configure dimension [%d] %s: " + e.getMessage(), i, names[i]);
            return null;
        }
    }
    // Check for a search
    boolean active = false;
    for (int i = 0; i < originalDimensions.length; i++) {
        if (originalDimensions[i].isActive()) {
            active = true;
            break;
        }
    }
    if (!active) {
        Utils.log(TITLE + " : No search range");
        return currentOptimum;
    }
    // Optionally use a reduced range (this is used for iteration)
    if (rangeReduction > 0 && rangeReduction < 1) {
        // Suppress dialogs and use the current settings
        nonInteractive = true;
        for (int i = 0; i < originalDimensions.length; i++) {
            double centre = point[i];
            double r = 0;
            if (originalDimensions[i].isActive()) {
                r = (originalDimensions[i].max - originalDimensions[i].min) * rangeReduction;
            }
            double lower = centre - r * 0.5;
            double upper = centre + r * 0.5;
            originalDimensions[i] = originalDimensions[i].create(lower, upper);
        }
    }
    analysisStopWatch = StopWatch.createStarted();
    // Store this for later debugging
    SearchResult<FilterScore> optimum = null;
    if (searchParam == 0 || searchParam == 2) {
        // Collect parameters for the range search algorithm
        pauseParameterTimer();
        boolean isStepSearch = searchParam == 2;
        // The step search should use a multi-dimension refinement and no range reduction
        SearchSpace.RefinementMode myRefinementMode = SearchSpace.RefinementMode.MULTI_DIMENSION;
        GenericDialog gd = null;
        boolean runAlgorithm = nonInteractive;
        if (!nonInteractive) {
            // Ask the user for the search parameters.
            gd = new GenericDialog(TITLE);
            gd.addMessage("Configure the " + SEARCH[searchParam] + " algorithm for " + ss_filter.getType());
            gd.addSlider("Width", 1, 5, pRangeSearchWidth);
            if (!isStepSearch) {
                gd.addNumericField("Max_iterations", pMaxIterations, 0);
                String[] modes = SettingsManager.getNames((Object[]) SearchSpace.RefinementMode.values());
                gd.addSlider("Reduce", 0.01, 0.99, pRangeSearchReduce);
                gd.addChoice("Refinement", modes, modes[pRefinementMode]);
            }
            gd.addNumericField("Seed_size", pSeedSize, 0);
            gd.showDialog();
            runAlgorithm = !gd.wasCanceled();
        }
        if (runAlgorithm) {
            SearchDimension[] dimensions = new SearchDimension[originalDimensions.length];
            if (!nonInteractive) {
                pRangeSearchWidth = (int) gd.getNextNumber();
                if (!isStepSearch) {
                    pMaxIterations = (int) gd.getNextNumber();
                    pRangeSearchReduce = gd.getNextNumber();
                    pRefinementMode = gd.getNextChoiceIndex();
                }
                pSeedSize = (int) gd.getNextNumber();
            }
            if (!isStepSearch)
                myRefinementMode = SearchSpace.RefinementMode.values()[pRefinementMode];
            for (int i = 0; i < dimensions.length; i++) {
                if (originalDimensions[i].isActive()) {
                    try {
                        dimensions[i] = originalDimensions[i].create(pRangeSearchWidth);
                        dimensions[i].setPad(true);
                        // Prevent range reduction so that the step search just does a single refinement step
                        dimensions[i].setReduceFactor((isStepSearch) ? 1 : pRangeSearchReduce);
                        // Centre on current optimum
                        dimensions[i].setCentre(point[i]);
                    } catch (IllegalArgumentException e) {
                        IJ.error(TITLE, String.format("Unable to configure dimension [%d] %s: " + e.getMessage(), i, names[i]));
                        return null;
                    }
                } else {
                    dimensions[i] = new SearchDimension(point[i]);
                }
            }
            // Check the number of combinations is OK
            long combinations = SearchSpace.countCombinations(dimensions);
            if (!nonInteractive && combinations > 10000) {
                gd = new GenericDialog(TITLE);
                gd.addMessage(String.format("%d combinations for the configured dimensions.\n \nClick 'Yes' to optimise.", combinations));
                gd.enableYesNoCancel();
                gd.hideCancelButton();
                gd.showDialog();
                if (!gd.wasOKed()) {
                    combinations = 0;
                }
            }
            if (combinations == 0) {
                resumeParameterTimer();
            } else {
                algorithm = SEARCH[searchParam] + " " + pRangeSearchWidth;
                ga_statusPrefix = algorithm + " " + ss_filter.getName() + " ... ";
                ga_iteration = 0;
                p_optimum = null;
                SearchSpace ss = new SearchSpace();
                ss.setTracker(this);
                if (pSeedSize > 0) {
                    // Add current optimum to seed
                    // Note: If we have an optimum and we are not seeding this should not matter as the dimensions 
                    // have been centred on the current optimum					
                    double[][] seed = new double[1][];
                    seed[0] = point;
                    // Sample without rounding as the seed will be rounded
                    double[][] sample = SearchSpace.sampleWithoutRounding(dimensions, pSeedSize - 1, null);
                    ss.seed(merge(sample, seed));
                }
                ConvergenceChecker<FilterScore> checker = new InterruptConvergenceChecker(0, 0, pMaxIterations);
                createGAWindow();
                resumeParameterTimer();
                optimum = ss.search(dimensions, new ParameterScoreFunction(), checker, myRefinementMode);
                if (optimum != null) {
                    // In case optimisation was stopped
                    IJ.resetEscape();
                    // Now update the parameters for final assessment
                    point = optimum.point;
                // Not required as the seed in now rounded
                //if (pSeedSize > 0)
                //{
                //	// The optimum may be off grid if it was from the seed
                //	point = enumerateMinInterval(point, names, originalDimensions);
                //}
                }
            }
        } else
            resumeParameterTimer();
    }
    if (searchParam == 1) {
        // Collect parameters for the enrichment search algorithm
        pauseParameterTimer();
        GenericDialog gd = null;
        boolean runAlgorithm = nonInteractive;
        if (!nonInteractive) {
            // Ask the user for the search parameters.
            gd = new GenericDialog(TITLE);
            gd.addMessage("Configure the " + SEARCH[searchParam] + " algorithm for " + ss_filter.getType());
            gd.addNumericField("Max_iterations", pMaxIterations, 0);
            gd.addNumericField("Converged_count", pConvergedCount, 0);
            gd.addNumericField("Samples", pEnrichmentSamples, 0);
            gd.addSlider("Fraction", 0.01, 0.99, pEnrichmentFraction);
            gd.addSlider("Padding", 0, 0.99, pEnrichmentPadding);
            gd.showDialog();
            runAlgorithm = !gd.wasCanceled();
        }
        if (runAlgorithm) {
            FixedDimension[] dimensions = Arrays.copyOf(originalDimensions, originalDimensions.length);
            if (!nonInteractive) {
                pMaxIterations = (int) gd.getNextNumber();
                pConvergedCount = (int) gd.getNextNumber();
                pEnrichmentSamples = (int) gd.getNextNumber();
                pEnrichmentFraction = gd.getNextNumber();
                pEnrichmentPadding = gd.getNextNumber();
            }
            algorithm = SEARCH[searchParam];
            ga_statusPrefix = algorithm + " " + ss_filter.getName() + " ... ";
            ga_iteration = 0;
            p_optimum = null;
            SearchSpace ss = new SearchSpace();
            ss.setTracker(this);
            // Add current optimum to seed
            double[][] seed = new double[1][];
            seed[0] = point;
            ss.seed(seed);
            ConvergenceChecker<FilterScore> checker = new InterruptConvergenceChecker(0, 0, pMaxIterations, pConvergedCount);
            createGAWindow();
            resumeParameterTimer();
            optimum = ss.enrichmentSearch(dimensions, new ParameterScoreFunction(), checker, pEnrichmentSamples, pEnrichmentFraction, pEnrichmentPadding);
            if (optimum != null) {
                // In case optimisation was stopped
                IJ.resetEscape();
                point = optimum.point;
            // Not required as the search now respects the min interval
            // Enumerate on the min interval to produce the final filter
            //point = enumerateMinInterval(point, names, originalDimensions);
            }
        } else
            resumeParameterTimer();
    }
    if (searchParam == 3) {
        // Collect parameters for the enumeration search algorithm
        pauseParameterTimer();
        SearchDimension[] dimensions = new SearchDimension[originalDimensions.length];
        for (int i = 0; i < dimensions.length; i++) {
            if (originalDimensions[i].isActive()) {
                try {
                    dimensions[i] = originalDimensions[i].create(0);
                } catch (IllegalArgumentException e) {
                    IJ.error(TITLE, String.format("Unable to configure dimension [%d] %s: " + e.getMessage(), i, names[i]));
                    return null;
                }
            } else {
                dimensions[i] = new SearchDimension(point[i]);
            }
        }
        GenericDialog gd = null;
        long combinations = SearchSpace.countCombinations(dimensions);
        if (!nonInteractive && combinations > 2000) {
            gd = new GenericDialog(TITLE);
            gd.addMessage(String.format("%d combinations for the configured dimensions.\n \nClick 'Yes' to optimise.", combinations));
            gd.enableYesNoCancel();
            gd.hideCancelButton();
            gd.showDialog();
            if (!gd.wasOKed()) {
                combinations = 0;
            }
        }
        if (combinations == 0) {
            resumeParameterTimer();
        } else {
            algorithm = SEARCH[searchParam];
            ga_statusPrefix = algorithm + " " + ss_filter.getName() + " ... ";
            ga_iteration = 0;
            p_optimum = null;
            SearchSpace ss = new SearchSpace();
            ss.setTracker(this);
            createGAWindow();
            resumeParameterTimer();
            optimum = ss.findOptimum(dimensions, new ParameterScoreFunction());
            if (optimum != null) {
                // In case optimisation was stopped
                IJ.resetEscape();
                // Now update the parameters for final assessment
                point = optimum.point;
            }
        }
    }
    IJ.showStatus("Analysing " + ss_filter.getName() + " ...");
    // Update the parameters using the optimum
    failCount = (int) Math.round(point[0]);
    residualsThreshold = sResidualsThreshold = point[1];
    duplicateDistance = point[2];
    // Refresh the coordinate store
    if (coordinateStore == null || duplicateDistance != coordinateStore.getResolution()) {
        coordinateStore = createCoordinateStore();
    }
    createResultsPrefix2();
    // (Re) Score the filter.
    // TODO - check this is now OK. Maybe remove the enumeration on the min interval grid
    // If scoring of filter here is different to scoring in the optimisation routine it is probably an ss_filter.clone() issue,
    // i.e. multi-threading use of the filter clone is not working.
    // Or it could be that the optimisation produced params off the min-interval grid
    FilterScoreResult scoreResult = scoreFilter(ss_filter);
    if (optimum != null) {
        if (scoreResult.score != optimum.score.score && scoreResult.criteria != optimum.score.criteria) {
            ParameterScoreResult r = scoreFilter((DirectFilter) ss_filter.clone(), minimalFilter, failCount, residualsThreshold, duplicateDistance, createCoordinateStore(duplicateDistance), false);
            System.out.printf("Weird re- score of the filter: %f!=%f or %f!=%f (%f:%f)\n", scoreResult.score, optimum.score.score, scoreResult.criteria, optimum.score.criteria, r.score, r.criteria);
        }
    }
    SimpleFilterScore max = new SimpleFilterScore(scoreResult, true, scoreResult.criteria >= minCriteria);
    analysisStopWatch.stop();
    if (showResultsTable) {
        BufferedTextWindow tw = null;
        if (resultsWindow != null)
            tw = new BufferedTextWindow(resultsWindow);
        addToResultsWindow(tw, scoreResult.text);
        if (resultsWindow != null)
            resultsWindow.getTextPanel().updateDisplay();
    }
    // Check the top result against the limits of the original dimensions
    StringBuilder sb = new StringBuilder(200);
    for (int j = 0; j < originalDimensions.length; j++) {
        if (!originalDimensions[j].isActive())
            continue;
        final double value = point[j];
        double lowerLimit = originalDimensions[j].getLower();
        double upperLimit = originalDimensions[j].getUpper();
        int c1 = Double.compare(value, lowerLimit);
        if (c1 <= 0) {
            sb.append(" : ").append(names[j]).append(' ').append(ComplexFilterScore.FLOOR).append('[').append(Utils.rounded(value));
            if (c1 == -1) {
                sb.append("<").append(Utils.rounded(lowerLimit));
            }
            sb.append("]");
        } else {
            int c2 = Double.compare(value, upperLimit);
            if (c2 >= 0) {
                sb.append(" : ").append(names[j]).append(' ').append(ComplexFilterScore.CEIL).append('[').append(Utils.rounded(value));
                if (c2 == 1) {
                    sb.append(">").append(Utils.rounded(upperLimit));
                }
                sb.append("]");
            }
        }
    }
    if (sb.length() > 0) {
        if (max.criteriaPassed) {
            Utils.log("Warning: Top filter (%s @ %s|%s) [%s] at the limit of the expanded range%s", ss_filter.getName(), Utils.rounded((invertScore) ? -max.score : max.score), Utils.rounded((invertCriteria) ? -minCriteria : minCriteria), limitFailCount + limitRange, sb.toString());
        } else {
            Utils.log("Warning: Top filter (%s @ -|%s) [%s] at the limit of the expanded range%s", ss_filter.getName(), Utils.rounded((invertCriteria) ? -max.criteria : max.criteria), limitFailCount + limitRange, sb.toString());
        }
    }
    // We may have no filters that pass the criteria
    String type = max.r.filter.getType();
    if (!max.criteriaPassed) {
        Utils.log("Warning: Filter does not pass the criteria: %s : Best = %s using %s", type, Utils.rounded((invertCriteria) ? -max.criteria : max.criteria), max.r.filter.getName());
        return null;
    }
    // Update without duplicates
    boolean allowDuplicates = false;
    // Re-use the atLimit and algorithm for the input optimum
    ComplexFilterScore newFilterScore = new ComplexFilterScore(max.r, currentOptimum.atLimit, currentOptimum.algorithm, currentOptimum.time, algorithm, analysisStopWatch.getTime());
    addBestFilter(type, allowDuplicates, newFilterScore);
    // Add spacer at end of each result set
    if (isHeadless) {
        if (showResultsTable)
            IJ.log("");
    } else {
        if (showResultsTable)
            resultsWindow.append("");
    }
    if (newFilterScore.compareTo(currentOptimum) <= 0)
        return newFilterScore;
    else {
        // Update the algorithm and time
        currentOptimum.paramAlgorithm = algorithm;
        currentOptimum.paramTime = analysisStopWatch.getTime();
    }
    return currentOptimum;
}
Also used : SearchSpace(gdsc.smlm.search.SearchSpace) GenericDialog(ij.gui.GenericDialog) NonBlockingGenericDialog(ij.gui.NonBlockingGenericDialog) BufferedTextWindow(gdsc.core.ij.BufferedTextWindow) SearchDimension(gdsc.smlm.search.SearchDimension) FixedDimension(gdsc.smlm.search.FixedDimension) FilterScore(gdsc.smlm.results.filter.FilterScore)

Aggregations

BufferedTextWindow (gdsc.core.ij.BufferedTextWindow)5 FilterScore (gdsc.smlm.results.filter.FilterScore)2 FixedDimension (gdsc.smlm.search.FixedDimension)2 SearchDimension (gdsc.smlm.search.SearchDimension)2 SearchSpace (gdsc.smlm.search.SearchSpace)2 GenericDialog (ij.gui.GenericDialog)2 NonBlockingGenericDialog (ij.gui.NonBlockingGenericDialog)2 TextWindow (ij.text.TextWindow)2 ArrayList (java.util.ArrayList)2 BasePoint (gdsc.core.match.BasePoint)1 FractionClassificationResult (gdsc.core.match.FractionClassificationResult)1 FastCorrelator (gdsc.core.utils.FastCorrelator)1 MaximaSpotFilter (gdsc.smlm.filters.MaximaSpotFilter)1 RampedSelectionStrategy (gdsc.smlm.ga.RampedSelectionStrategy)1 SimpleMutator (gdsc.smlm.ga.SimpleMutator)1 SimpleRecombiner (gdsc.smlm.ga.SimpleRecombiner)1 SimpleSelectionStrategy (gdsc.smlm.ga.SimpleSelectionStrategy)1 PeakResultPoint (gdsc.smlm.ij.plugins.ResultsMatchCalculator.PeakResultPoint)1 DirectFilter (gdsc.smlm.results.filter.DirectFilter)1 Filter (gdsc.smlm.results.filter.Filter)1