FindClusters for 1D numeric lists: ten methods, three count modes - #54
Open
msollami wants to merge 3 commits into
Open
FindClusters for 1D numeric lists: ten methods, three count modes#54msollami wants to merge 3 commits into
msollami wants to merge 3 commits into
Conversation
Adds the two-argument Nearest: the element(s) of list at minimum
Abs[element - x], as a List in original order.
All tied elements are returned, which is what fixes the algorithm.
Nearest[{1, 5, 10}, 3] is {1, 5}, not {1}. The shape is MinimalBy's
(src/sort.c:663-716) — find the minimum, then collect every distance
equal to it — rather than the RankedMin quickselect, whose comparator
carries an original-index tiebreak (src/sort.c:932) that exists to make
ties impossible and would return a single element. Input order among
ties falls out of the ascending collect pass, so there is no tie logic
in the file.
Distance composes the existing internal_subtract and internal_abs the
way comparisons.c:313-314 already does; no distance helper was added.
Abs of a complex difference is its modulus, so Nearest[{3 + 4 I, 1}, 0]
is {1} with no extra code.
Distances are ordered by numeric value, not by expr_compare. Canonical
order is wrong here in both directions, and each direction breaks the
all-ties guarantee. It settles a value tie between different ExprTypes
on the type enum (src/sort.c:376), so Nearest[{0, 2.0}, 1] would answer
{0} — distances 1 and 1.0 are equal but Integer sorts before Real —
dropping a tied element from the one function whose contract is to
return them all. And for atoms that are not both integer-like it
compares get_numeric_value() doubles (src/sort.c:372-377), so
Nearest[{1/3, 1/3 + 1/10^18}, 0] would report a tie between distances
that differ exactly. Subtracting is exact where it must be and inexact
only where the input already was: 1 - 1.0 is 0.0, a genuine tie, while
1/3 - (1/3 + 1/10^18) is the exact Rational[-1, 10^18]. nearest_sign
reports undecidable separately from zero, which expr_numeric_sign
cannot — it returns a bare 0 for both and recognises neither MPFR nor a
bigint-component Rational.
Nearest diverges from MinimalBy in one respect, deliberately: every
distance must be a real number or the call stays unevaluated.
MinimalBy[{1, a, 3}, Abs[# - 2] &] answers {1, 3}, dropping the symbolic
element because expr_compare orders symbols after all numbers — a
plausible wrong answer. Gating on the distance rather than the element
covers a symbolic element, a symbolic target and a non-real complex in
one check.
Two limitations are pinned by test rows rather than left silent. A
symbolic real such as Pi declines instead of being numericalized. A
rational with a bigint component declines for a reason outside this
file: builtin_abs does not evaluate one (Abs[1/1000] is 1/1000,
Abs[1/10^25] is Abs[1/10^25], and Sign has the same gap), so the
distance arrives unevaluated and the gate rejects it. Both rows flip the
day the underlying behaviour changes.
Only the two-argument form lands here; the n-nearest, radius, rule,
all-pairs and NearestTo forms and the DistanceFunction option are
separate. A packed list is materialised on the way in since Nearest is
not on pack.c's AWARE list; a visible NDArray is not a List and stays
unevaluated rather than being silently truncated.
New: src/list/nearest.{c,h}, registered in list_init.c with Protected.
30 acceptance rows in tests/test_list.c::test_nearest. make check-c99,
check-packed-aware and check-array-exactness pass; leak-free under a
differential leaks run (0 bytes at 200 and 20000 iterations across the
success, gate-bail, mixed-type-tie and decline paths).
The distance comparator was intransitive, and an intransitive comparator
feeding a sort produces order-dependent output.
Two distances of differing exactness were ordered by subtracting them,
but internal_subtract widens the pair to a double and loses the exact
operand. With plain int64 values, no bigints or rationals needed:
a = 2^60, b = 2.0^60, c = 2^60 + 1
a - b -> 0.0 so the comparator said a == b
c - b -> 0.0 so it said c == b
a - c -> -1 so it said a < c
all three at once. Nearest[{2^60 + 1, 2.0^60}, 2^60] therefore reported
a tie between an element at distance 1 and one at distance 0.
A double is exactly a rational, m * 2^e, so lifting both sides into mpq
and comparing there is exact and gives a total order. nearest_to_mpq is
shaped after mod_quot_expr_to_mpq (core.c:2569), which is static there
and has no Real case; the Real case is the point of the new helper.
Non-finite Reals and MPFR operands are not liftable and fall through to
the subtraction path, which is what handled them before.
This changes one user-visible behaviour, deliberately, and Mathematica
agrees on both halves of it: 1 and 1.0 still tie, because 1.0 lifts to
exactly 1, while 0.1 and 1/10 no longer do. Mathematica's
Nearest[{0.9, 11/10}, 1] is {0.9}, not both, and now so is ours.
Nearest[{0, 2.0}, 1] and Nearest[{1.5, 5/2}, 2] still return both
elements, since those reals are exactly representable.
Found by code review, not by the table: all 29 existing rows passed
throughout, because none of them probed a mixed-exactness comparison
near the limits of double precision. Four rows added that do.
make check-c99 passes; 33 rows in test_nearest pass; leak-free at 20000
iterations across the success, tie, exact-rational and decline paths.
Partitions a 1D numeric list into clusters, as a list of lists. Clusters
come out in order of the first occurrence of a member and elements keep
their input order, matching Gather's convention.
THREE COUNT MODES, NOT TWO. n and UpTo[n] are genuinely different, and
Mathematica gives them separate error classes with different capability
lists. n forces exactly that many clusters; UpTo[n] returns the natural
count when it is already at or below n. On {1,2,10,12,3,1,13,25}, whose
natural count is three, n -> 4 splits {10,12,13} while UpTo[4] does not.
Both cap at the distinct-value count, since no method can separate two
equal elements.
The 10x3 capability matrix is a static table checked before dispatch:
KMeans and KMedoids need a count, the five density methods need
Automatic, and Spectral takes Automatic and UpTo[n] but not a bare n.
That last cell is transcribed from the runtime's own allowed-lists,
which contradict the documentation bullets -- the bullets place Spectral
in neither constraint list, implying it takes a fixed count, and the
runtime rejects it.
WE DO NOT REPRODUCE MATHEMATICA'S OUTPUT, and cannot. It auto-selects a
distance function and preprocesses the data by unpublished rules, so
even single-linkage with an explicit count disagrees with the textbook
algorithm: on {1,4,9,16,25,36}, gaps 3,5,7,9,11, it cuts 11 and 7 rather
than the two largest. Its Automatic count is not a gap rule either.
Each method here is the textbook algorithm with stated semantics, and
the acceptance table is the specification. Three divergences are pinned
as rows rather than left implicit.
Agglomerate and SpanningTree share one implementation because in 1D they
are the same computation: single-linkage equals cutting the widest edges
of the minimum spanning tree, and on a line the spanning tree is the
sorted adjacency chain. A row requires the two names to agree.
EXACT WHERE IT MATTERS. Ordering, distinctness and fixed-count gap
selection all compare the elements themselves through list_numeric_cmp,
never a double projection, so FindClusters[{1/10^25, 2/10^25, 1}, 2]
works -- which Nearest cannot do, because it routes its distance through
Abs and builtin_abs declines on a bigint-component rational
(complex.c:418-421). A sorted 1D pass never needs Abs: the gap between
sorted neighbours is non-negative by construction. A row fails if anyone
rewrites a gap as Abs[b - a].
Supporting change: the numeric comparator written for Nearest was
file-static and FindClusters is its second caller, so it moved to
list_common.{c,h} as list_real_number_q / list_numeric_sign /
list_numeric_cmp. The 33 Nearest rows pass untouched.
Findings from code review, all fixed here and each covered by a row:
- fc_emit_clusters' OOM unwind double-freed the failing slot. expr_copy
is a refcount bump, so the second release would recycle a node the
caller still owned.
- JarvisPatrick returned four clusters for nine copies of 7, and DBSCAN
split a duplicated pair when MinPoints exceeded its multiplicity. The
equal-elements invariant is now enforced once, centrally, on exact
zero gaps, so no method can reintroduce it.
- KMeans/KMedoids returned fewer than a fixed n on tie-heavy data, and
on data whose distinctness is invisible to doubles: {2^60, 2^60+1,
2^60+2} with n=3 gave one cluster. Seeding and boundary selection are
now exact.
- The Automatic threshold averaged the two middle gaps, which makes a
split arithmetically impossible for two or three elements --
{0, 1, 10^12} was one cluster. It now takes the lower median,
computed on the gap Exprs, which also removes the double projection
and the scale cliff where 10^25 worked and 10^400 did not.
- MeanShift fragmented evenly spaced data (Range[40] gave 34 clusters)
because the merge tolerance was tied to the bandwidth; the two scales
are now separate. UpTo[n] was identical to a bare n for KMeans and
KMedoids. The two quadratic methods now decline above 4000 rather
than appearing to hang.
110 acceptance rows, every expected value produced by the built binary
rather than predicted, plus three invariants asserted as loops:
partition (no element lost or duplicated), exactly-n, and
equal-elements. make check-c99, check-packed-aware and
check-array-exactness pass; the fastpath sweep records FindClusters in
OFF_BUFFER rather than half-building a buffer path; leak-free at 200 and
20 000 iterations across every method and every decline path.
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
1 participant
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
FindClusters[list],[list, n]and[list, UpTo[n]]partition a 1D numeric list into clusters, as a list of lists. All ten of Mathematica's namedMethodvalues are implemented.Clusters come out in order of the first occurrence of a member; elements keep their input order.
Three count modes, not two
nandUpTo[n]are genuinely different — Mathematica gives them separate error classes with different capability lists:The 10×3 capability matrix is a static table checked before dispatch:
AutomaticUpTo[n]n"Agglomerate","SpanningTree""KMeans","KMedoids""Spectral""DBSCAN","GaussianMixture","JarvisPatrick","MeanShift","NeighborhoodContraction"The
Spectralrow is transcribed from the runtime's own allowed-lists, which contradict the documentation bullets — the bullets placeSpectralin neither constraint list, implying it takes a fixed count, and the runtime rejects it. Thirty rows cover the matrix, one per cell.We do not reproduce Mathematica's output, and cannot
It auto-selects a distance function and preprocesses the data by unpublished rules — its own messages name the choice — so even single-linkage with an explicit count disagrees with the textbook algorithm:
gaps
3,5,7,9,11, and it cut11and7rather than the two largest. ItsAutomaticcount is not a gap rule either. So each method here is the textbook algorithm with stated semantics, and the acceptance table is the specification. No success criterion is "matches Mathematica"; three divergences are pinned as rows.AgglomerateandSpanningTreeshare one implementation because in 1D they are the same computation — single-linkage equals cutting the widest MST edges, and on a line the spanning tree is the sorted adjacency chain. A row requires the two names to agree.Exact where it matters
Ordering, distinctness and fixed-count gap selection all compare the elements themselves, never a double projection:
Nearestcannot do this — it routes its distance throughAbs, andbuiltin_absdeclines on a bigint-component rational (complex.c:418-421). A sorted 1D pass never needsAbs: the gap between sorted neighbours is non-negative by construction. A row fails if anyone rewrites a gap asAbs[b - a].Supporting change
The numeric comparator written for
Nearestwas file-static and this is its second caller, so it moved tolist_common.{c,h}aslist_real_number_q/list_numeric_sign/list_numeric_cmp. The 33Nearestrows pass untouched.Fixed during review
Each is covered by a row; none was reachable from the original table, which is why review found them and the table did not.
fc_emit_clustersdouble-free on the OOM unwind.expr_copyis a refcount bump, so the second release would recycle a node the caller still owned. The only memory-safety defect.JarvisPatrickreturned four clusters for nine copies of7;DBSCANsplit a duplicated pair whenMinPointsexceeded its multiplicity. Now enforced once, centrally, on exact zero gaps, so no method can reintroduce it.ncame back short on tie-heavy data and where distinctness is invisible to doubles:{2^60, 2^60+1, 2^60+2}withn=3gave one cluster. Seeding and boundary selection are now exact.Automaticcould not split 2–3 element lists. Averaging the two middle gaps makes it arithmetically impossible — the largest gap is one of the averaged values.{0, 1, 10^12}was one cluster. Now the lower median, computed on the gapExprs, which also removes the scale cliff where10^25worked and10^400did not.MeanShiftfragmented evenly spaced data (Range[40]→ 34 clusters) because the merge tolerance was tied to the bandwidth; the two scales are now separate.UpTo[n]was identical to a barenforKMeans/KMedoids.MeanShiftwas ~16 min at 10⁵).Testing
110 acceptance rows — every expected value produced by the built binary, none predicted — plus three invariants asserted as loops: partition (no element lost or duplicated), exactly-
n, and equal-elements.make check-c99,check-packed-aware,check-array-exactnesspass. The fastpath sweep recordsFindClustersinOFF_BUFFERrather than half-building a buffer path. Leak-free at 200 and 20 000 iterations across every method and every decline path.Known limitations, deliberately deferred
Automaticmedian uses an insertion sort — correct but O(n²); wantsnd_select_kth. Default path is ~0.65 s at 10⁵.MeanShift/NeighborhoodContractionare quadratic; capped rather than reworked onto the sorted-window primitives.CriterionFunctionandPerformanceGoalare accepted and have no effect — stated in the docs rather than left implicit.