Example 1 with FixedDimension
use of gdsc.smlm.search.FixedDimension 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 FixedDimension
use of gdsc.smlm.search.FixedDimension 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)2
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
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
DirectFilter (gdsc.smlm.results.filter.DirectFilter)1
Filter (gdsc.smlm.results.filter.Filter)1
FilterSet (gdsc.smlm.results.filter.FilterSet)1
IDirectFilter (gdsc.smlm.results.filter.IDirectFilter)1
MultiPathFilter (gdsc.smlm.results.filter.MultiPathFilter)1
ParameterType (gdsc.smlm.results.filter.ParameterType)1
ArrayList (java.util.ArrayList)1
RandomDataGenerator (org.apache.commons.math3.random.RandomDataGenerator)1
