fix(fractal-dim): oracle-correct box-count & perimeter dimensions#366
fix(fractal-dim): oracle-correct box-count & perimeter dimensions#366sameeul wants to merge 8 commits into
Conversation
Smallest power of two >= a, unlike closest_pow2 which is strictly greater. Used to pad an ROI to a box-counting grid without wasting an octave (ceil_pow2(256)=256 vs 512). Co-Authored-By: Claude Opus 4.8 (1M context) <noreply@anthropic.com>
Box counting is grid-registration sensitive. Power2PaddedImageMatrix (and its oversized _NT variant) padded to a power of two strictly larger than the ROI and centered it, so the ROI never lined up with the coarse boxes: halving the box size did not cleanly quadruple the count and the dimension came out systematically low (a filled square read 1.75 instead of 2.0, Sierpinski 1.39 instead of 1.585). Pad tight (ceil_pow2) and place the ROI at the grid origin. Validated against the ImageJ/FracLac box-count oracle. Co-Authored-By: Claude Opus 4.8 (1M context) <noreply@anthropic.com>
Box counting (calculate_boxcount_fdim): count boxes directly from the ROI pixel cloud on a tight, origin-aligned grid, and auto-switch on ROI size - large ROIs use a single grid via the padded mask matrix (fast early-exit tile scan), small ROIs shift the grid over 4 origins and take the minimum count (as FracLac does), which removes the residual registration bias that only bites when few box sizes fit. Recovers analytic shapes exactly (square->2, line->1, Sierpinski->log2(3)) and reproduces the ImageJ/FracLac shifting-grid oracle on real ROIs. Estimator: replace the old mean-of-local-slopes (an endpoint slope) with loglog_slope, a least-squares fit of log(measure) vs log(scale) shared by both features (box-count D = -slope, perimeter D = 1 - slope). Perimeter (calculate_perimeter_fdim): rewrite as a clean closed-contour Richardson divider - perimeter as a double (was truncated to int), drop the broken tail chord (it added a variable-length closing segment and used squared distance instead of the Euclidean chord), and scale by the actual mean ruler length. Recovers the Koch snowflake (log4/log3=1.2619) and a smooth disk (1.0); the single-vertex walk matches analytic truth (a multi-start average was tried and rejected - it degraded the Koch recovery to 1.37). Co-Authored-By: Claude Opus 4.8 (1M context) <noreply@anthropic.com>
…I tests Validate both fractal features against external oracles: - Analytic shapes (tests/python/test_fractal_dim_oracle.py): box-count square->2, line->1, Sierpinski->log2(3); perimeter disk->1, Koch snowflake->log4/log3. numpy-only, runs in CI. - Arbitrary 154-px ROI: box-count vs the ImageJ shifting-grid oracle (same method, 1.389, tight) and the divider perimeter vs box-count-of-edge (cross-method, 1.163, agree ~6%). Move the C++ shape2d fractal assertions to the 154-px pixelIntensityFeaturesTestData ROI (test_shape2d_fractal_dimension_154px_oracle): it fits 5 box sizes, enough for BOTH features to be oracle-validated - box-count same-method (1% tol) and perimeter cross-method (10% tol). The 26-px morphology fixture only fits 3 scales, where the perimeter has no convergent cross-method oracle, so the fractal features are no longer asserted there. Other shape2d goldens (convex hull, extrema, circles, ...) stay on the 26-px fixture. Co-Authored-By: Claude Opus 4.8 (1M context) <noreply@anthropic.com>
a946076 to
38f6989
Compare
| // and biased the dimension low (a filled square read 1.75 instead of 2.0). | ||
| int bigSide = std::max(aabb.get_width(), aabb.get_height()); | ||
| StatsInt paddedSide = Nyxus::closest_pow2 (bigSide); | ||
| StatsInt paddedSide = Nyxus::ceil_pow2 (bigSide); |
There was a problem hiding this comment.
Nyxus::closest_pow2 is not used anymore, can be removed
| inline int ceil_pow2(const int a) | ||
| { | ||
| int p = 1; | ||
| while (p < a) |
There was a problem hiding this comment.
Overflow issue for a > 2^30, p <<= 1 goes negative and goes into infinite loop
Suggesting either rewrite it like this
inline uint32_t ceil_pow2(uint32_t a) {
if (a == 0) return 1;
if (a > (1u << 31)) return 0; // or UINT_MAX / assert — pick a policy
uint32_t x = a - 1;
x |= x>>1; x|=x>>2; x|=x>>4; x|=x>>8; x|=x>>16;
return x + 1;
}
or switch to uint64_t
#include
inline uint64_t ceil_pow2(uint32_t a) {
return a <= 1 ? 1 : std::bit_ceil<uint64_t>(a);
}
| std::vector<std::pair<double, double>> coverage; | ||
|
|
||
| coverage.push_back({ s, cnt }); | ||
| if (paddedSide > 32) |
There was a problem hiding this comment.
Magic number, need a comment on why we look for 32
| // boundary dimension); Nyxus' Richardson divider agrees to ~6% (cross-method | ||
| // convergent validity), asserted with a 10% tolerance -> 1.163 | ||
| static std::unordered_map<std::string, double> oracle_fractal_154px_golden_values{ | ||
| {"FRACT_DIM_BOXCOUNT", 1.3891}, |
There was a problem hiding this comment.
Nit (non-blocking): flag that 1.3891 is a same-method pin, not external ground truth.
This value comes from the ImageJ/FracLac shifting-grid box count — the same algorithm Nyxus runs — so the 1% assertion below is really a regression/self-consistency pin, not an independent oracle. That's fine and the header comment (L72) already says "SAME method as Nyxus", but since it sits in a map literally named ...oracle... next to the genuinely cross-method FRACT_DIM_PERIMETER, it's easy to misread 1.3891 as external corroboration. Suggest a one-line inline marker so the distinction survives a skim, e.g.:
{"FRACT_DIM_BOXCOUNT", 1.3891}, // self-consistency pin: SAME shifting-grid method as Nyxus (not independent validation)
{"FRACT_DIM_PERIMETER", 1.163}, // cross-method (edge box-count vs divider) — the real external check
The independent validation here is the analytic shapes (square→2, line→1, Sierpinski→log₂3, …) and the cross-method perimeter; the 154-px box-count is a pin.
There was a problem hiding this comment.
Updated the test with a new image of arbitrary ROI, typically found in microscopy images. The numbers have been vetted by running ImageJ/FracLac locally.
| # algorithm than nyxus' Richardson divider, but estimates the SAME boundary-complexity | ||
| # dimension (box dim == compass/divider dim for boundary curves). The two agree to | ||
| # 0.06 here - cross-method convergent validity -> 1.163. | ||
| ORACLE_BOXCOUNT_154 = 1.389 # same-method oracle; nyxus computes 1.389 (<0.001) |
There was a problem hiding this comment.
Nit (non-blocking): same as the C++ side — mark this as a same-method pin.
ORACLE_BOXCOUNT_154 = 1.389 is the offline ImageJ shifting-grid box count, i.e. the same method Nyxus uses, so abs(bc - 1.389) < 0.02 is a self-consistency check rather than independent corroboration (the trailing comment "<0.001" makes that clear — a near-exact match is expected precisely because it's the same algorithm). Consider renaming to make intent obvious at the call site, or at minimum keep the comment explicit:
ORACLE_BOXCOUNT_154 = 1.389 # SAME-method pin (ImageJ shifting-grid == nyxus); not independent
ORACLE_PERIMETER_154_EDGE = 1.163 # cross-method (edge box-count vs divider) — independent check
No change to the assertion itself needed; the tight 0.02 tolerance is correct for a same-method pin. The genuinely independent Python checks remain the analytic-shape tests and the cross-method perimeter.
There was a problem hiding this comment.
Updated the test with a new image of arbitrary ROI, typically found in microscopy images. The numbers have been vetted by running ImageJ/FracLac locally.
|
|
||
| std::vector<Pixel2> K; | ||
| r.merge_multicontour(K); | ||
|
|
There was a problem hiding this comment.
It is worth checking that K is truly single countour
// The divider assumes ONE closed contour. Multi-contour ROIs (holes / disconnected
// fragments) concatenate into K and inject spurious seam chords that bias the dimension.
assert(r.contour_count() <= 1 && "perimeter fdim assumes single contour");
| // Fractal dimensions are validated on the irregular 154-px pixelIntensityFeaturesTestData ROI, which | ||
| // fits enough box sizes for both features to be oracle-validated (the 26-px morphology fixture only | ||
| // fits 3 scales, where the perimeter has no convergent cross-method oracle). | ||
| static void calculate_fractal_154px_feature_values(std::vector<std::vector<double>>& fvals) |
There was a problem hiding this comment.
Nit (non-blocking, maintenance watch-point): this harness hand-copies the production ROI-preprocessing path.
load_test_roi_data only fills raw_pixels, so L339–342 re-implement inline what the real masked loader does — make_nonanisotropic_aabb() → aux_image_matrix.allocate(...) → calculate_from_pixelcloud(...) → initialize_fvals() (plus the hardcoded label 102). That's correct for the current loader and the oracle values are meaningful, but nothing keeps it in sync: if production ever changes how an ROI is finalized (AABB/anisotropy convention, matrix fill, a new normalization step), this fixture keeps building the ROI the old way and would still pass — validating fractal dimensions against an ROI the real pipeline would never produce.
Two ways to de-risk, in order of preference:
route the test through the same finalize entry point the loader uses instead of duplicating its steps, so a production change flows in automatically; or
if no shared entry point exists, add a // KEEP IN SYNC WITH marker next to the make_nonanisotropic_aabb block and an assert pinning the resulting AABB/matrix dims, so drift is visible. Not a blocker for this PR (the reconstruction is already commented and matches production today) — just flagging the hidden coupling so whoever touches ROI preprocessing later updates it here too.
Replace the inline 154-px ROI fractal oracle with a large (512x512) irregular ROI committed as tests/data/fractal_blob512_seg.ome.tif, loaded from disk by both the C++ and Python suites. Box-count is validated same-method (FracLac single origin-aligned grid, 1.8706, 1% tol); perimeter cross-method (FracLac edge box-count, 1.0493, 10% tol). Co-Authored-By: Claude Opus 4.8 (1M context) <noreply@anthropic.com>
- helpers.h: rewrite ceil_pow2 branch-free — the `while (p < a) p <<= 1` shift loop overflowed to a negative int (infinite loop) for a > 2^30; remove now-unused closest_pow2. - image_matrix.h: drop the stale closest_pow2 reference in the padding comment. - fractal_dim.cpp: document the box-count large/small-ROI threshold (32) and the single-contour assumption of the divider perimeter. - tests: mark the box-count golden as a same-method self-consistency pin (vs the cross-method perimeter check) in both suites, and flag the hand-copied ROI-finalize block as KEEP IN SYNC with the production loader. Co-Authored-By: Claude Opus 4.8 (1M context) <noreply@anthropic.com>
…cy pin The box-count oracle (1.8706) is ImageJ/FracLac's output — an independent implementation of the same box-count method, so Nyxus' value is vetted by a third-party tool (it caught the grid-registration bug a self-pin never would). Same method, independent implementation; the perimeter check is cross-method. Co-Authored-By: Claude Opus 4.8 (1M context) <noreply@anthropic.com>
…er helper Comment-7 follow-up: instead of hand-copying the ROI finalize (AABB + image matrix) inline with a KEEP IN SYNC marker, use load_masked_test_roi_data - the shared helper the other 2D shape tests already use - so the finalize path lives in one place. Co-Authored-By: Claude Opus 4.8 (1M context) <noreply@anthropic.com>
Summary
Corrects the accuracy of the 2D fractal-dimension features (
FRACT_DIM_BOXCOUNT,FRACT_DIM_PERIMETER) and validates them against external oracles — analytic shapes with a knowndimension, and the ImageJ/FracLac tool. Follows up the sign/aggregation fix already on
main.Before → after (the accuracy bugs, on shapes with known truth):
Root cause
Box counting is grid-registration sensitive.
Power2PaddedImageMatrixpadded to a power of twostrictly larger than the ROI (
closest_pow2) and centered it, so the ROI never lined upwith the coarse boxes — halving the box size didn't cleanly quadruple the count, biasing the
dimension low on every shape. The perimeter path additionally truncated the perimeter to
intand added a broken tail chord (variable length + squared distance instead of the Euclidean chord),
pushing
FRACT_DIM_PERIMETERout of the valid[1,2]range.Changes
Feature (
fractal_dim.cpp,image_matrix.h,helpers/helpers.h)ceil_pow2) with the ROI at the grid origin (replacesclosest_pow2+ centering). Auto-switch on ROI size — large ROIs use a single aligned grid (fastearly-exit tile scan), small ROIs shift the grid over 4 origins and take the minimum count (as
FracLac does). Perf-neutral (large-ROI ~16 ms, small ~35 µs).
doubleperimeter, dropped tail chord,mean-ruler-length scale.
loglog_slopeestimator (box-slope, perimeter1-slope).ceil_pow2rewritten branch-free / overflow-safe; the now-unusedclosest_pow2removed.Tests (
tests/python/test_fractal_dim_oracle.py,tests/test_shape_morphology_2d.h,tests/data/)Sierpinski→log2(3); perimeter disk→1, Koch→log4/log3.
tests/data/fractal_blob512_seg.ome.tif— a large (512×512) irregularsingle ROI, the kind of object a microscopy segmentation produces — loaded from disk by both
the C++ and Python suites and validated against offline ImageJ/FracLac oracles.
Oracle validation
Each row is a real validation, not a self-referential snapshot:
algorithm. This validates Nyxus' implementation (a third-party tool agreeing is genuine
corroboration — it is what caught the 1.75-vs-2.0 bug; a self-consistency pin never would). Same
method ⇒ tight tolerance expected.
validity (~10%).
Why a large ROI. A fractal dimension is a scaling measurement — meaningful only across enough
scales. The 512² blob fits 9 dyadic box sizes and exercises the production single-grid (large-ROI)
path. On it the grid registration is stable: Nyxus' single grid and FracLac's origin-aligned grid
agree to ~0.5% — a genuine implementation check. (An earlier iteration validated on a 154-px ROI,
which fit only ~5 box sizes and sat in the small-ROI shifting-grid path — a weaker, borderline
fixture.) FracLac is GUI-only, so the box-count oracle (single origin-aligned grid) and the perimeter
oracle (box count of the ROI edge) were read from the FracLac GUI offline and pinned as constants.
Research note: divider vs box-counting for
FRACT_DIM_PERIMETERFRACT_DIM_PERIMETERuses the Richardson divider; the perimeter oracle is FracLac's box-count of theROI edge — a different algorithm. They estimate the same boundary dimension:
divider/compass, Hausdorff dimensions all equal) — which is why Nyxus' divider recovers the Koch
snowflake to 0.004. Divergence appears only on finite, digitized, non-self-similar curves and
shrinks with resolution (on the 512² blob, 9 box sizes, divider and box-count-of-edge agree to ~1%).
box-count 1.143 vs divider 1.130 (Δ≈0.013) and states "the divider dimension more accurately
represents the irregularity of a coastline"
(Sci. Rep. 2021).
method); box-counting was adopted mainly because it removes the divider's bookkeeping ambiguities
and automates for any set — not because it is more accurate for lines. See Klinkenberg,
A review of methods used to determine the fractal dimension of linear features
(Math. Geology 1994).
We keep the divider (its DIN ISO 9276-6 "structured walk" heritage) because it recovers the one
analytic ground truth essentially exactly. A multi-start averaging variant (the divider analog of
shifting grids) was implemented and rejected — it degraded the Koch recovery from 1.266 to 1.37.
🤖 Generated with Claude Code