Add support for polydisperse features

[putnokiabel]

Goal

Currently, trackpy only allows setting a single diameter (monodisperse features) for detection. This works if all features are the same size, but does not work if an image/video has features of multiple sizes (polydisperse) which is still a common setup that is currently not supported.
The workarounds at the moment are:

• track everything with the diameter of the largest feature (works only if the size differences between features are minimal and the distance between features is large), or
• track the image multiple times with different diameters and combine the results by deduplicating features that are close to each other (takes significantly more compute time, deduplication is not trivial, accuracy is mediocre).

This PR introduces true support for polydisperse features by allowing the user to set diameter=Polydisperse(min_diameter, max_diameter) and making changes to tracking to facilitate this range of diameters at every step (preprocessing, refining, deduplication, spiff, etc).
The existing tracking for monodisperse particles (i.e. if you do not set diameter=Polydisperse) remains completely unchanged.

Benchmark:

I added a benchmark (run using python benchmarks/polydisperse_grid.py) that tests a combination of different scenarios:

• size differences (3x, 5x, 10x)
• feature density (low, med, high)
• noise level (low, med, high)

and compares the baseline (tracking everything with the max diameter) vs polydisperse (spiff enabled for both).
We compare:

• Recall r (the proportion of features that get identified at all, within a tolerance of 1.5px distance)
• Precision p (the proportion of features found that are real features)
• Error e (RMS of distance between actual and estimated feature positions, for the features that were identified correctly)
• Total number of correct features identified
• Total tracking time

Results

2D benchmark

range dens  noise    N | baseline             poly                
----------------------------------------------------------------------------
x3    low   low    125 | r1.00 e0.020 p1.00  r1.00 e0.038 p1.00
x3    low   med    120 | r1.00 e0.045 p1.00  r1.00 e0.052 p1.00
x3    low   high   121 | r1.00 e0.089 p1.00  r1.00 e0.088 p1.00
x3    med   low    300 | r1.00 e0.022 p1.00  r1.00 e0.040 p1.00
x3    med   med    300 | r1.00 e0.035 p1.00  r1.00 e0.051 p1.00
x3    med   high   300 | r1.00 e0.083 p1.00  r1.00 e0.087 p1.00
x3    high  low    300 | r0.82 e0.020 p1.00  r1.00 e0.040 p1.00
x3    high  med    300 | r0.82 e0.036 p1.00  r1.00 e0.050 p1.00
x3    high  high   300 | r0.83 e0.075 p1.00  r1.00 e0.085 p1.00
x5    low   low    102 | r1.00 e0.018 p1.00  r1.00 e0.033 p1.00
x5    low   med    105 | r1.00 e0.057 p1.00  r1.00 e0.060 p1.00
x5    low   high   113 | r1.00 e0.162 p0.90  r1.00 e0.087 p0.97
x5    med   low    270 | r0.54 e0.027 p1.00  r1.00 e0.047 p1.00
x5    med   med    266 | r0.59 e0.039 p1.00  r1.00 e0.058 p1.00
x5    med   high   266 | r0.59 e0.113 p0.99  r1.00 e0.097 p1.00
x5    high  low    300 | r0.44 e0.065 p1.00  r1.00 e0.071 p1.00
x5    high  med    300 | r0.47 e0.074 p0.99  r1.00 e0.072 p1.00
x5    high  high   300 | r0.40 e0.097 p0.97  r1.00 e0.103 p1.00
x10   low   low     82 | r0.49 e0.009 p1.00  r1.00 e0.037 p1.00
x10   low   med     91 | r0.46 e0.046 p1.00  r1.00 e0.058 p1.00
x10   low   high    87 | r0.45 e0.238 p1.00  r1.00 e0.093 p1.00
x10   med   low    186 | r0.18 e0.080 p0.97  r1.00 e0.071 p1.00
x10   med   med    186 | r0.17 e0.195 p0.86  r1.00 e0.080 p1.00
x10   med   high   186 | r0.20 e0.216 p0.90  r1.00 e0.113 p1.00
x10   high  low    300 | r0.07 e0.327 p0.67  r1.00 e0.110 p1.00
x10   high  med    300 | r0.08 e0.349 p0.70  r1.00 e0.105 p1.00
x10   high  high   300 | r0.04 e0.304 p0.44  r1.00 e0.123 p1.00
----------------------------------------------------------------------------
Total correct (within 1.5px) of 5906 true:  baseline=3463  poly=5906  (poly/baseline=1.71x)
Total grid time:  baseline=0.58s  poly=0.92s  (poly/baseline=1.59x)

As you can see, baseline does OK in the 3x range when the features aren't too dense, but fails completely in other scenarios (recall drops dramatically), whereas Polydisperse has (near) perfect recall and precision across every simulated scenario.
Note that Polydisperse has slightly higher tracking errors than baseline, however it's still the same order of magnitude and still gives good subpixel accuracy.

Polydisperse is, on a per-feature-found basis, actually faster than baseline (it finds 71% more features but only takes 59% longer).

3D anisotropic benchmark (mimics confocal microscopy)

confocal aspect=(0.3333333333333333, 1, 1)  shape=(48, 256, 256)  noise_std=12
range dens     N | baseline             poly                
----------------------------------------------------------------------
x5    low    150 | r0.89 e0.458 p0.88  r0.99 e0.300 p1.00
x5    med    150 | r0.66 e0.221 p0.79  r1.00 e0.291 p1.00
x5    high   150 | r0.63 e0.250 p0.73  r1.00 e0.306 p1.00
x10   low     76 | r0.21 e0.348 p0.89  r1.00 e0.412 p0.99
x10   med    150 | r0.09 e0.442 p1.00  r0.99 e0.313 p1.00
x10   high   150 | r0.08 e0.165 p0.92  r0.99 e0.308 p1.00
----------------------------------------------------------------------
Total correct (within 1.5px) of 826 true:  baseline=370  poly=822  (poly/baseline=2.22x)
Total detect+refine time:  baseline=2.00s  poly=5.63s  (poly/baseline=2.81x)

In this benchmark you again see that poly has (near) perfect recall and precision across the grid, while baseline fails.
Poly identifies 2.22x as many features and takes 2.81x as long, so is slightly slower than baseline (but is still a great trade-off).

Current scope:

• Support in locate, batch and helper methods (like refine and spiff). This PR does not include linking (which could be significantly improved if we use particle size for linking, assuming features don't change in size over time).
• Support for both 2d and 3d.
• Support for both isotropic and anisotropic images.
• Support up to 10x range in feature size.

Future work:

• Use particle size in linking (see also Assigning sizes to each particle. #694)
• Add multiprocessing (some parts of the updating tracking process could benefit from additional parallellization/multiprocessing, as some steps are ran independently for different batches of diameters)
• Test for scenarios where there is more than a 10x range in feature sizes, and diagnose bottlenecks
• Test the accuracy of the output diameter. If needed, do a second pass of (re)determining the diameters after refining COM and SPIFF, for better diameter accuracy.
• Add function to show distribution in feature diameters of a certain image or video. This could also help users determine the right tracking diameter even in monodisperse images.
• Test if recalculating feature diameter in-between COM refine iterations would improve accuracy.

[putnokiabel]

Note:
Currently I only enabled this in polydisperse mode. However, I think this could lead to a performance improvement also in the existing, monodisperse setup.
By doing this piece of filtering upfront, we can skip wasteful refine cycles for multiple peaks that are touching (and thus clearly part of the same feature).

[putnokiabel]

Note:
This is not only a nice-to-have output for users of the package, it's also required to make SPIFF work with multiple feature sizes.

[b-grimaud]

Very excited for this, we're tracking biological objects of varying sizes and setting a proper diameter is challenging. We miss a lot of meaningful signal if the diameter is too small, but linking (and MSD-based filtering) times become impractical if we set it too high.

[putnokiabel]

Very excited for this, we're tracking biological objects of varying sizes and setting a proper diameter is challenging. We miss a lot of meaningful signal if the diameter is too small, but linking (and MSD-based filtering) times become impractical if we set it too high.

@b-grimaud
That's great! Feel free to test this PR on your dataset and leave feedback (since I've only tested using synthetic data so far) if you find any gaps :)

[putnokiabel]

Update: I added support for anisotropic images (useful for confocal images, which is the usecase that prompted me to work in this issue) so the scope is a bit broader now.

@nkeim curious what you think!

[nkeim]

This is great!! Thanks for this PR. I'm sorry that I've been traveling this whole month but I will try to take a closer look in the next couple weeks.

I want to put out a release before the end of the summer—we already have some big improvements! Depending on circumstances and how much testing folks can give it, this may stay in the development branch until the next release. My group will be eager to test it out!