Skip to content

Keep floating point error in normalized weights from adding a replicate in weighted_quantile(n = ...) - #278

Open
haomeng797-ship-it wants to merge 1 commit into
mjskay:masterfrom
haomeng797-ship-it:fix-n-rep-rounding
Open

Keep floating point error in normalized weights from adding a replicate in weighted_quantile(n = ...)#278
haomeng797-ship-it wants to merge 1 commit into
mjskay:masterfrom
haomeng797-ship-it:fix-n-rep-rounding

Conversation

@haomeng797-ship-it

Copy link
Copy Markdown

Found while working through #267. This is a separate, narrower bug, so I've split it out rather than mixing it into that discussion.

The problem

The docs for n promise that summarising a sample with duplicates into a weighted sample and then passing the original sample size gives back quantile() on the original sample. That breaks on some inputs:

x = rep(1:3, times = c(10, 8, 7))

quantile(x, 0.75, type = 4, names = FALSE)
#> [1] 2.75
weighted_quantile(1:3, 0.75, weights = c(10, 8, 7), n = 25, type = 4, names = FALSE)
#> [1] 2.857143

The weights are normalized to weights / sum(weights), then expanded with n_rep = ceiling(weights * n). For the third group that round trip is not exact:

7/25 * 25
#> [1] 7
print(7/25 * 25, digits = 17)
#> [1] 7.0000000000000009

so ceiling() returns 8 instead of 7. That group gets an extra replicate, the total becomes 26 rather than 25, and every p_k downstream shifts.

Sweeping 500 random integer samples across types 4–9 (3000 comparisons), 1.0% disagree with quantile() on the original sample. The existing test for this correspondence happens to use counts 4:1, whose normalization is exact, so it doesn't catch it.

The fix

Shrink by a relative tolerance before rounding up, so weights that are exact multiples of 1/n stay put:

n_rep = pmax(1, ceiling(weights * n * (1 - sqrt(.Machine$double.eps))))

pmax(1, ...) keeps very small weights from being shrunk to zero replicates, which would drop the point entirely. After the change the same 3000-comparison sweep has no mismatches.

Tests

Added a regression test alongside the existing equivalence test, using counts that do trigger the rounding. It fails on main for types 4, 5, 7, and 9, and passes with the fix. test.weighted_quantile.R, test.weighted_ecdf.R, test.weighted_hist.R, and test.point_interval.R all pass.

Worth noting the tolerance is a pragmatic fix at the point of rounding, not a change to how weights are represented. If you'd rather have the expansion work from unnormalized weights (where counts stay exact), that would be a more structural fix, and I'm happy to redo it that way.

weighted_quantile(n = ...) expands each point into ceiling(weights * n)
replicates. When the weights come from counts, normalizing and multiplying
back by n can land a few ulp above a whole number, so a point picks up an
extra replicate and the result no longer matches quantile() on the original
sample, which is what n is documented to guarantee.

Found while looking at mjskay#267.
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

None yet

Projects

None yet

Development

Successfully merging this pull request may close these issues.

1 participant