Repository navigation
Rework micrometry quad matching - #72
Merged
Merged
Conversation
Micrometry could fail without saying so on two-camera bead data: - Mirrored channels. The quad labelling rules (xc <= xd, xc + xd <= 1) pick one labelling for each group of four points. A mirror image often breaks them, so the swapped-code comparison in isMatch() found only 21% of mirrored quads with 200 beads, and none with 50. The 'other' quads are now also made from the positions with x negated, and fit using the original positions, so the affine transform includes the mirroring. - Different magnifications. min_size and max_size were applied in each camera's own pixels, so at a scale of 1.25 only 42% of the quads existed in both cameras. Quads now come from each localization and its 8 nearest neighbors, with no size limits by default. min_size and max_size remain as optional limits on |AB|. - The score. log(fg_p / bg_p), with bg_p = 1 / (movie pixels) and a fixed 1 px width, meant something different for each movie size and localization count. The score is now -ln of the chance that a wrong transform pairs this many localizations one to one within r = max(1, 3 sqrt(2) localization_error): a Poisson tail, multiplied by the number of transforms tried. A score above 10 means the same thing for any field. Pairing one to one also stops several localizations scoring against the same reference. - Tolerance. A box tolerance of 0.01 is too tight for 0.5 to 1 px of localization error. Candidates are now checked nearest code first, up to a budget of 2000, with a tolerance of 0.04. - Refinement. The transform returned and saved was the fit to the four points of one quad. It is now refined by least squares on all the one to one pairs until they stop changing. - Speed. Quads are built with numpy, matched through a KD-tree on their codes, and the candidate transforms are fit and ranked in one batch. With 200 beads, building the quads takes 0.03 s rather than 4 s, and the whole match 0.25 s; comparing the quads alone took about 16 s. With 2000 beads on 2048 x 2048 a match takes 3.6 s. (One process on a slow core.) - Micrometry() no longer fails when not given max_neighbors. findTransformFromPoints() takes two (N, 2) arrays rather than HDF5 files, and saveMapping() writes the map file that runMicrometry() saves, in the same format. runMicrometry(), Micrometry.findTransform() and the command line keep their arguments, with new defaults (max_neighbors 8, tolerance 0.04, no size limits) and two new ones, localization_error (default 0.5 px) and candidate_budget. Measured on 9,600 simulated trials: 50 or 200 beads on 512 x 512, camera 2 rotated up to 10 degrees, scaled 0.8 to 1.25, shifted up to 30 px and mirrored in half the trials, each bead missed with probability 0 to 0.7 in each camera, 0.1 to 1 px localization error: - 88% of transforms are within max(1, 3 sqrt(2) sigma) px RMS across the field, the same as an independent implementation of this pipeline (87.7%). On the 480 trials a port of the previous code was run on, 78% against 50%. - Mirroring makes no difference (87.9% against 88.0%), nor does the camera scale, and with up to 30% of beads missed 1 px of error does as well as 0.1 px. - On 5,160 pairs of unrelated point lists (10 to 2000 points, 256 to 2048 px fields, 0.1 to 1 px error) the highest score was 7.3. Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
Micrometry reads SA HDF5 localizations, but make_test_data.py wrote Insight3 files and the refine_map.py command line read them, so neither worked with it. makeTestData() writes two HDF5 files, and can now rotate, scale, shift and mirror the second field and add localization error. Localizations moved outside the field are dropped. refineMap() refines an existing map using the same one to one pairing and least squares as micrometry, with the pairing radius set by localization_error rather than a fixed 1 px. It is for updating a known map, for example with beads imaged on another day, as micrometry now refines the maps it finds. The new tests run runMicrometry() from HDF5 files on mirrored, rotated and scaled data, check the saved map and that the HDF5 and array entry points agree, and refine a map that is 1.5 px off. Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
analysis.rst suggested tuning --max_size for sparse beads, which no longer applies by default. It now describes the score and when to raise --localization_error. The multicolor diagnostic uses runMicrometry()'s defaults. On its simulated beads the three maps score 50 to 55, where they scored 11.0 against the threshold of 10 before, and are within 0.05 px RMS of the simulated shifts, as before. Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
The score says that the pairs are not chance, not that the map is accurate everywhere. With few paired beads, or beads covering only part of the field, the map can be several pixels off away from them. In the simulated trials, 20 of 8,235 accepted maps were outside the success radius across the field. In 16 of them almost every pair was a real bead pair, but the affine fit to beads covering 10 to 50% of the field did not extrapolate. The other 4 had chance pairs pulling the fit, at 1 px localization error with heavy dropout. mappingError() estimates the RMS error of the map across the field from the least squares fit: the residual scatter of the pairs, times the leverage of the fit averaged over the reference field. On the accepted simulated maps the actual error is a median 1.0x this estimate (90th percentile 1.5x, 99th 2.2x). Micrometry keeps it (getMappingError(), getPairCount()) and prints it when verbose. runMicrometry() warns when the estimate is more than a third of the pairing radius, which is the error of one pair of localizations unless the radius is at its 1 px minimum. It still saves the map. In the simulated trials this warns about 15 of the 20 inaccurate maps and 165 of the 8,215 accurate ones (2.0%), the latter all with heavy dropout and a median error of 43% of the radius. Warning on fewer than 10 pairs as well catches no more and doubles the false warnings. On the multicolor diagnostic's beads the estimate is 0.01 to 0.04 px, in line with the measured 0.01 to 0.05 px, and there is no warning. Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
Add this suggestion to a batch that can be applied as a single commit.This suggestion is invalid because no changes were made to the code.Suggestions cannot be applied while the pull request is closed.Suggestions cannot be applied while viewing a subset of changes.Only one suggestion per line can be applied in a batch.Add this suggestion to a batch that can be applied as a single commit.Applying suggestions on deleted lines is not supported.You must change the existing code in this line in order to create a valid suggestion.Outdated suggestions cannot be applied.This suggestion has been applied or marked resolved.Suggestions cannot be applied from pending reviews.Suggestions cannot be applied on multi-line comments.Suggestions cannot be applied while the pull request is queued to merge.Suggestion cannot be applied right now. Please check back later.
Summary
Micrometry now finds the map between two bead localization lists when the channels are mirrored, differ in magnification, or have 0.5 to 1 px of localization error. It saves a map refined on all the paired beads, and reports a score that means the same thing for any field size and bead count.
localization_error), and the transform is refit by least squares until the pairs stop changing.runMicrometry()warns when the map's estimated error across the field is more than a third of the pairing radius (see below). The map is still saved.findTransformFromPoints()takes two (N, 2) arrays, andsaveMapping()writes the map file.refine_mapandmake_test_datanow read and write HDF5.runMicrometry(),Micrometry.findTransform(), the command line and the map file format keep their arguments. The defaults change (max_neighbors8,tolerance0.04, no size limits), and there are two new arguments,localization_error(0.5 px) andcandidate_budget(2000).quads.MicroQuadandmicrometry.fgProbability()are gone.Results
Simulated two-camera trials (9,600): 50 or 200 beads on 512 × 512. Camera 2 is rotated up to 10°, scaled 0.8 to 1.25, shifted up to 30 px, and mirrored in half the trials. Each bead is missed with probability p = 0 to 0.7 in each camera, with σ = 0.1 to 1 px of localization error. A trial succeeds if the transform is within max(1, 3√2σ) px RMS of the true one across the field.
Success does not change with camera scale.
Unrelated point lists (5,160 pairs; 10 to 2000 points, 256 to 2048 px fields, σ = 0.1 to 1): the highest score was 7.3.
Time for one match from HDF5 files, mirrored data at scale 1.1, one process on a slow core:
The two versions' scores are on different scales; for both, 10 is the acceptance threshold.
Multicolor diagnostic: the three maps score 50 to 55, where they scored 11.0 before, and are within 0.05 px RMS of the simulated shifts, as before.
Accurate pairs, inaccurate map
The score tests whether the pairs could have arisen by chance, not how well the map extrapolates. In 20 of the 8,235 accepted simulation trials, the transform was outside the success radius across the field. All 20 have heavy dropout (p ≥ 0.3), and only 7 to 18 beads were paired.
The last commit adds
mappingError(), which estimates the map's RMS error across the field from the least squares fit: the pairs' residual scatter times the fit's leverage, averaged over the reference field. On accepted maps the actual error is a median 1.0× the estimate (1.5× at the 90th percentile, 2.2× at the 99th).Micrometrykeeps it (getMappingError(),getPairCount()) and prints it when verbose.runMicrometry()warns when the estimate is more than a third of the pairing radius, which equals the error of a single pair of localizations unless the radius is at its 1 px minimum. On the accepted simulation maps:The accurate maps it warns about all have heavy dropout, and their median error is 43% of the radius. On the multicolor diagnostic's beads the estimate is 0.01 to 0.04 px, in line with the measured 0.01 to 0.05 px, and there is no warning.
Testing
pytest: 288 passed (Linux, Python 3.14).test_micrometry.pygrows from 2 tests to 8: quad codes against the scalar definition; mirror + scale 1.2 + σ = 0.5 + 30% dropout; unrelated lists rejected across field sizes and counts;runMicrometry()end to end from HDF5; refining a map; the error estimate and warning for beads spread over the field and clustered in a corner.🤖 Generated with Claude Code