Examples with LocalList uk.ac.sussex.gdsc.core.utils.LocalList used on opensource projects

Example 41 with LocalList

use of uk.ac.sussex.gdsc.core.utils.LocalList in project GDSC-SMLM by aherbert.

the class TcPalmAnalysis method createCumulativeCountData.

/**
 * Creates the cumulative count data.
 *
 * @param createPlotData set to true to create the plot data arrays
 * @return the cumulative count data
 */
private CumulativeCountData createCumulativeCountData(LocalList<ClusterData> clusters, boolean createPlotData) {
    final TIntIntHashMap all = new TIntIntHashMap(maxT - minT + 1);
    clusters.forEach(c -> c.results.forEach(peak -> all.adjustOrPutValue(peak.getFrame(), 1, 1)));
    final int[] frames = all.keys();
    final int[] counts = all.values();
    SortUtils.sortData(counts, frames, true, false);
    return new CumulativeCountData(frames, counts, createPlotData);
}

Example 42 with LocalList

use of uk.ac.sussex.gdsc.core.utils.LocalList in project GDSC-SMLM by aherbert.

the class TcPalmAnalysis method reportAnalysis.

/**
 * Report statistics on the analysis results.
 *
 * @param settings the settings
 * @param clusters the clusters
 * @param calibration the data calibration
 */
private static void reportAnalysis(TcPalmAnalysisSettings settings, LocalList<ClusterData> clusters, DataCalibration calibration) {
    final WindowOrganiser wo = new WindowOrganiser();
    final Consumer<HistogramPlotBuilder> action = builder -> {
    /* noop. */
    };
    plotHistogram(settings.getShowSizeHistogram(), wo, clusters, "Size", c -> c.results.size(), null, action);
    plotHistogram(settings.getShowDurationHistogram(), wo, clusters, "Duration (" + calibration.getTimeUnitName() + ")", c -> calibration.timeConverter.convert(c.getDuration()), null, action);
    plotHistogram(settings.getShowAreaHistogram(), wo, clusters, "Area (" + calibration.getDistanceUnitName() + "^2)", c -> calibration.convertArea(c.getArea()), null, action);
    plotHistogram(settings.getShowDensityHistogram(), wo, clusters, "Density (" + calibration.getDistanceUnitName() + "^-2)", c -> c.results.size() / calibration.convertArea(c.getArea()), Double::isFinite, action);
    wo.tile();
}

Example 43 with LocalList

use of uk.ac.sussex.gdsc.core.utils.LocalList in project GDSC-SMLM by aherbert.

the class TraceDiffusion method filterTraces.

/**
 * Split traces to contiguous traces and filter traces that are not the minimum length. Re-assigns
 * the ID for the output traces.
 *
 * @param name the name
 * @param traces the traces
 * @param minimumTraceLength the minimum trace length
 * @param ignoreEnds the ignore ends
 * @return The new traces
 */
private static Trace[] filterTraces(String name, Trace[] traces, int minimumTraceLength, boolean ignoreEnds) {
    final LocalList<Trace> list = new LocalList<>(traces.length);
    final int minLength = (ignoreEnds) ? minimumTraceLength + 2 : minimumTraceLength;
    final Consumer<Trace> action = (ignoreEnds) ? t -> {
        if (t.size() >= minLength) {
            t.removeEnds();
            list.add(t);
            t.setId(list.size());
        }
    } : t -> {
        if (t.size() >= minLength) {
            list.add(t);
            t.setId(list.size());
        }
    };
    final Consumer<ArrayPeakResultStore> action2 = (ignoreEnds) ? r -> {
        if (r.size() >= minLength) {
            r.remove(0);
            r.remove(r.size() - 1);
            final Trace t = new Trace(r);
            list.add(t);
            t.setId(list.size());
        }
    } : r -> {
        if (r.size() >= minLength) {
            final Trace t = new Trace(r);
            list.add(t);
            t.setId(list.size());
        }
    };
    final ArrayPeakResultStore results = new ArrayPeakResultStore(11);
    for (final Trace trace : traces) {
        if (trace.size() < minLength) {
            // Too short
            continue;
        }
        if (trace.size() == trace.getTail().getFrame() - trace.getHead().getFrame() + 1) {
            // Contiguous
            action.accept(trace);
        } else {
            // Split the trace
            int t1 = trace.getHead().getFrame();
            for (int i = 0; i < trace.size(); i++) {
                final PeakResult peak = trace.get(i);
                final int t2 = peak.getFrame();
                if (t2 - t1 > 1) {
                    // Non-contiguous
                    action2.accept(results);
                    results.clear();
                }
                t1 = t2;
                results.add(peak);
            }
            // Final trace
            action2.accept(results);
            results.clear();
        }
    }
    ImageJUtils.log("Filtered results '%s' : %s split and filtered to %d using " + "minimum length %d (Ignore ends = %b)", name, TextUtils.pleural(traces.length, "trace"), list.size(), minimumTraceLength, ignoreEnds);
    return list.toArray(new Trace[0]);
}

Example 44 with LocalList

use of uk.ac.sussex.gdsc.core.utils.LocalList in project GDSC-SMLM by aherbert.

the class TrackPopulationAnalysis method run.

@Override
public void run(String arg) {
    SmlmUsageTracker.recordPlugin(this.getClass(), arg);
    if (MemoryPeakResults.isMemoryEmpty()) {
        IJ.error(TITLE, "No localisations in memory");
        return;
    }
    settings = Settings.load();
    // Saved by reference so just save now
    settings.save();
    // Read in multiple traced datasets
    // All datasets must have the same pixel pitch and exposure time
    // Get parameters
    // Convert datasets to tracks
    // For each track compute the 4 local track features using the configured window
    // 
    // Optional:
    // Fit a multi-variate Gaussian mixture model to the data
    // (using the configured number of components/populations)
    // Assign each point in the track using the model.
    // Smooth the assignments.
    // 
    // The alternative is to use the localisation category to assign populations.
    // 
    // Plot histograms of each track parameter, coloured by component
    final List<MemoryPeakResults> combinedResults = new LocalList<>();
    if (!showInputDialog(combinedResults)) {
        return;
    }
    final boolean hasCategory = showHasCategoryDialog(combinedResults);
    if (!showDialog(hasCategory)) {
        return;
    }
    ImageJUtils.log(TITLE + "...");
    final List<Trace> tracks = getTracks(combinedResults, settings.window, settings.minTrackLength);
    if (tracks.isEmpty()) {
        IJ.error(TITLE, "No tracks. Please check the input data and min track length setting.");
        return;
    }
    final Calibration cal = combinedResults.get(0).getCalibration();
    final CalibrationReader cr = new CalibrationReader(cal);
    // Use micrometer / second
    final TypeConverter<DistanceUnit> distanceConverter = cr.getDistanceConverter(DistanceUnit.UM);
    final double exposureTime = cr.getExposureTime() / 1000.0;
    final Pair<int[], double[][]> trackData = extractTrackData(tracks, distanceConverter, exposureTime, hasCategory);
    final double[][] data = trackData.getValue();
    // Histogram the raw data.
    final Array2DRowRealMatrix raw = new Array2DRowRealMatrix(data, false);
    final WindowOrganiser wo = new WindowOrganiser();
    // Store the histogram data for plotting the components
    final double[][] columns = new double[FEATURE_NAMES.length][];
    final double[][] limits = new double[FEATURE_NAMES.length][];
    // Get column data
    for (int i = 0; i < FEATURE_NAMES.length; i++) {
        columns[i] = raw.getColumn(i);
        if (i == FEATURE_D) {
            // Plot using a logarithmic scale
            SimpleArrayUtils.apply(columns[i], Math::log10);
        }
        limits[i] = MathUtils.limits(columns[i]);
    }
    // Compute histogram bins
    final int[] bins = new int[FEATURE_NAMES.length];
    if (settings.histogramBins > 0) {
        Arrays.fill(bins, settings.histogramBins);
    } else {
        for (int i = 0; i < FEATURE_NAMES.length; i++) {
            bins[i] = HistogramPlot.getBins(StoredData.create(columns[i]), BinMethod.FD);
        }
        // Use the maximum so all histograms look the same
        Arrays.fill(bins, MathUtils.max(bins));
    }
    // Compute plots
    final Plot[] plots = new Plot[FEATURE_NAMES.length];
    for (int i = 0; i < FEATURE_NAMES.length; i++) {
        final double[][] hist = HistogramPlot.calcHistogram(columns[i], limits[i][0], limits[i][1], bins[i]);
        plots[i] = new Plot(TITLE + " " + FEATURE_NAMES[i], getFeatureLabel(i, i == FEATURE_D), "Frequency");
        plots[i].addPoints(hist[0], hist[1], Plot.BAR);
        ImageJUtils.display(plots[i].getTitle(), plots[i], ImageJUtils.NO_TO_FRONT, wo);
    }
    wo.tile();
    // The component for each data point
    int[] component;
    // The number of components
    int numComponents;
    // Data used to fit the Gaussian mixture model
    double[][] fitData;
    // The fitted model
    MixtureMultivariateGaussianDistribution model;
    if (hasCategory) {
        // Use the category as the component.
        // No fit data and no output model
        fitData = null;
        model = null;
        // The component is stored at the end of the raw track data.
        final int end = data[0].length - 1;
        component = Arrays.stream(data).mapToInt(d -> (int) d[end]).toArray();
        numComponents = MathUtils.max(component) + 1;
        // In the EM algorithm the probability of each data point is computed and normalised to
        // sum to 1. The normalised probabilities are averaged to create the weights.
        // Note the probability of each data point uses the previous weight and the algorithm
        // iterates.
        // This is not a fitted model but the input model so use
        // zero weights to indicate no fitting was performed.
        final double[] weights = new double[numComponents];
        // Remove the trailing component to show the 'model' in a table.
        createModelTable(Arrays.stream(data).map(d -> Arrays.copyOf(d, end)).toArray(double[][]::new), weights, component);
    } else {
        // Multivariate Gaussian mixture EM
        // Provide option to not use the anomalous exponent in the population mix.
        int sortDimension = SORT_DIMENSION;
        if (settings.ignoreAlpha) {
            // Remove index 0. This shifts the sort dimension.
            sortDimension--;
            fitData = Arrays.stream(data).map(d -> Arrays.copyOfRange(d, 1, d.length)).toArray(double[][]::new);
        } else {
            fitData = SimpleArrayUtils.deepCopy(data);
        }
        final MultivariateGaussianMixtureExpectationMaximization mixed = fitGaussianMixture(fitData, sortDimension);
        if (mixed == null) {
            IJ.error(TITLE, "Failed to fit a mixture model");
            return;
        }
        model = sortComponents(mixed.getFittedModel(), sortDimension);
        // For the best model, assign to the most likely population.
        component = assignData(fitData, model);
        // Table of the final model using the original data (i.e. not normalised)
        final double[] weights = model.getWeights();
        numComponents = weights.length;
        createModelTable(data, weights, component);
    }
    // Output coloured histograms of the populations.
    final LUT lut = LutHelper.createLut(settings.lutIndex);
    IntFunction<Color> colourMap;
    if (LutHelper.getColour(lut, 0).equals(Color.BLACK)) {
        colourMap = i -> LutHelper.getNonZeroColour(lut, i, 0, numComponents - 1);
    } else {
        colourMap = i -> LutHelper.getColour(lut, i, 0, numComponents - 1);
    }
    for (int i = 0; i < FEATURE_NAMES.length; i++) {
        // Extract the data for each component
        final double[] col = columns[i];
        final Plot plot = plots[i];
        for (int n = 0; n < numComponents; n++) {
            final StoredData feature = new StoredData();
            for (int j = 0; j < component.length; j++) {
                if (component[j] == n) {
                    feature.add(col[j]);
                }
            }
            if (feature.size() == 0) {
                continue;
            }
            final double[][] hist = HistogramPlot.calcHistogram(feature.values(), limits[i][0], limits[i][1], bins[i]);
            // Colour the points
            plot.setColor(colourMap.apply(n));
            plot.addPoints(hist[0], hist[1], Plot.BAR);
        }
        plot.updateImage();
    }
    createTrackDataTable(tracks, trackData, fitData, model, component, cal, colourMap);
// Analysis.
// Assign the original localisations to their track component.
// Q. What about the start/end not covered by the window?
// Save tracks as a dataset labelled with the sub-track ID?
// Output for the bound component and free components track parameters.
// Compute dwell times.
// Other ...
// Track analysis plugin:
// Extract all continuous segments of the same component.
// Produce MSD plot with error bars.
// Fit using FBM model.
}

Example 45 with LocalList

use of uk.ac.sussex.gdsc.core.utils.LocalList in project GDSC-SMLM by aherbert.

the class MultivariateGaussianMixtureExpectationMaximization method createMixed.

/**
 * Helper method to create a multivariate Gaussian mixture model which can be used to initialize
 * {@link #fit(MixtureMultivariateGaussianDistribution)}.
 *
 * <p>This method can be used with the data supplied to the instance constructor to try to
 * determine a good mixture model at which to start the fit, but it is not guaranteed to supply a
 * model which will find the optimal solution or even converge.
 *
 * <p>The weights for each component will be uniform.
 *
 * @param data data for the distribution
 * @param component the component for each data point
 * @return Multivariate Gaussian mixture model estimated from the data
 * @throws IllegalArgumentException if data has less than 2 rows, if there is a size mismatch
 *         between the data and components length, or if the number of components is less than 2.
 */
public static MixtureMultivariateGaussianDistribution createMixed(double[][] data, int[] component) {
    ValidationUtils.checkArgument(data.length >= 2, "Estimation requires at least 2 data points: %d", data.length);
    ValidationUtils.checkArgument(data.length == component.length, "Data and component size mismatch: %d != %d", data.length, component.length);
    final int numRows = data.length;
    final int numCols = data[0].length;
    ValidationUtils.checkArgument(numCols >= 2, "Multivariate Gaussian requires at least 2 data columns: %d", numCols);
    // Sort the data
    final ClassifiedData[] sortedData = new ClassifiedData[data.length];
    for (int i = 0; i < numRows; i++) {
        sortedData[i] = new ClassifiedData(data[i], component[i]);
    }
    Arrays.sort(sortedData, (o1, o2) -> Integer.compare(o1.value, o2.value));
    ValidationUtils.checkArgument(sortedData[0].value != sortedData[sortedData.length - 1].value, "Mixture model requires at least 2 data components");
    // components of mixture model to be created
    final LocalList<MultivariateGaussianDistribution> distributions = new LocalList<>();
    int from = 0;
    while (from < sortedData.length) {
        // Find the end
        int to = from + 1;
        final int comp = sortedData[from].value;
        while (to < sortedData.length && sortedData[to].value == comp) {
            to++;
        }
        // number of data records that will be in this component
        final int numCompRows = to - from;
        // data for this component
        final double[][] compData = new double[numCompRows][];
        // mean of each column for the data
        final double[] columnMeans = new double[numCols];
        // populate and create component
        int count = 0;
        for (int i = from; i < to; i++) {
            final double[] values = sortedData[i].data;
            compData[count++] = values;
            for (int j = 0; j < numCols; j++) {
                columnMeans[j] += values[j];
            }
        }
        SimpleArrayUtils.multiply(columnMeans, 1.0 / numCompRows);
        distributions.add(new MultivariateGaussianDistribution(columnMeans, covariance(columnMeans, compData)));
        from = to;
    }
    // uniform weight for each bin
    final int numComponents = distributions.size();
    final double[] weights = SimpleArrayUtils.newDoubleArray(numComponents, 1.0 / numComponents);
    return new MixtureMultivariateGaussianDistribution(weights, distributions.toArray(new MultivariateGaussianDistribution[0]));
}