use of gdsc.smlm.search.SearchDimension 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;
}
use of gdsc.smlm.search.SearchDimension in project GDSC-SMLM by aherbert.
the class BenchmarkFilterAnalysis method enumerateMinInterval.
/**
* Enumerate on the min interval to convert an off grid result to one on the grid
*
* @param best
* The optimum
* @param enabled
* Array specifying which parameters are enabled
* @param indices
* Array specifying which parameter indices to search
* @return The optimum on the min interval grid
*/
private Chromosome<FilterScore> enumerateMinInterval(Chromosome<FilterScore> best, boolean[] enabled, int[] indices) {
// Enumerate on the min interval to produce the final filter
ss_filter = (DirectFilter) best;
es_optimum = null;
SearchDimension[] dimensions2 = new SearchDimension[ss_filter.getNumberOfParameters()];
for (int i = 0; i < indices.length; i++) {
int j = indices[i];
if (enabled[j]) {
double minIncrement = ss_filter.getParameterIncrement(j);
try {
double value = ss_filter.getParameterValue(j);
double max = Maths.ceil(value, minIncrement);
double min = Maths.floor(value, minIncrement);
dimensions2[i] = new SearchDimension(min, max, minIncrement, 1);
dimensions2[i].setCentre(value);
dimensions2[i].setIncrement(minIncrement);
} catch (IllegalArgumentException e) {
IJ.error(TITLE, String.format("Unable to configure dimension [%d] %s: " + e.getMessage(), j, ss_filter.getParameterName(j)));
dimensions2 = null;
break;
}
}
}
if (dimensions2 != null) {
// Add dimensions that have been missed
for (int i = 0; i < dimensions2.length; i++) if (dimensions2[i] == null)
dimensions2[i] = new SearchDimension(ss_filter.getParameterValue(i));
SearchSpace ss = new SearchSpace();
ss.setTracker(this);
SearchResult<FilterScore> optimum = ss.findOptimum(dimensions2, this);
if (optimum != null) {
best = ((SimpleFilterScore) optimum.score).r.filter;
}
}
return best;
}
use of gdsc.smlm.search.SearchDimension in project GDSC-SMLM by aherbert.
the class BenchmarkFilterAnalysis method enumerateMinInterval.
@SuppressWarnings("unused")
private double[] enumerateMinInterval(double[] point, String[] names, FixedDimension[] originalDimensions) {
p_optimum = null;
SearchDimension[] dimensions2 = new SearchDimension[originalDimensions.length];
for (int i = 0; i < dimensions2.length; i++) {
if (originalDimensions[i].isActive()) {
double minIncrement = originalDimensions[i].minIncrement;
try {
double value = point[i];
double max = Maths.ceil(value, minIncrement);
double min = Maths.floor(value, minIncrement);
dimensions2[i] = new SearchDimension(min, max, minIncrement, 1);
dimensions2[i].setCentre(value);
dimensions2[i].setIncrement(minIncrement);
} catch (IllegalArgumentException e) {
IJ.error(TITLE, String.format("Unable to configure dimension [%d] %s: " + e.getMessage(), i, names[i]));
dimensions2 = null;
break;
}
} else
dimensions2[i] = new SearchDimension(point[i]);
}
if (dimensions2 != null) {
SearchSpace ss = new SearchSpace();
ss.setTracker(this);
SearchResult<FilterScore> optimum = ss.findOptimum(dimensions2, new ParameterScoreFunction());
if (optimum != null) {
point = optimum.point;
}
}
return point;
}
use of gdsc.smlm.search.SearchDimension 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;
}
Aggregations