Skip to content

fix: use log1p for small p in stats/base/dists/planck/quantile - #15764

Open
Abhist17 wants to merge 2 commits into
stdlib-js:developfrom
Abhist17:fix/planck-quantile-small-p
Open

Abhist17 wants to merge 2 commits into
stdlib-js:developfrom
Abhist17:fix/planck-quantile-small-p

Conversation

@Abhist17

@Abhist17 Abhist17 commented Oct 1, 2026 •

Copy link
Copy Markdown
Contributor

Description

What is the purpose of this pull request?

This pull request:

  • uses -log1p( -p ) instead of -ln( 1.0 - p ) in stats/base/dists/planck/quantile (JS, factory and C).
  • clamps the result at 0, like geometric/quantile, so p = 0 and inputs where -log1p( -p ) / lambda underflows (e.g. quantile( 5.0e-324, 2.0 )) return 0 instead of -1.

Forming 1 - p in double precision rounds away the low-order bits of p. That gives two kinds of wrong answer:

var quantile = require( '@stdlib/stats/base/dists/planck/quantile' );

// `1 - p` rounds to `1`, so the result is `ceil(0) - 1`, outside the support:
quantile( 1.0e-20, 1.0 );
// develop: -1
// this PR: 0

// `1 - 1e-15` is off by about 11%, which a small `lambda` turns into whole steps:
quantile( 1.0e-15, 1.0e-17 );
// develop: 111
// this PR: 100 (2048-bit evaluation of `ceil(-ln(1-p)/lambda) - 1` agrees)

SciPy's planck.ppf already uses log1p(-q), so the existing Python fixtures pass unchanged with exact equality.

Tests: a new "small p" test in test.quantile.js, test.factory.js and test.native.js covers p = 1e-20 and 1e-300 (expected 0) and two small-lambda cases (expected 100, checked against 2048-bit BigFloat). Four of those assertions fail on develop.

Related Issues

Does this pull request have any related issues?

None. Same bug as #15762 (weibull/quantile) and #15763 (exponential/quantile).

Questions

Any questions for reviewers of this pull request?

No.

Other

Any other information relevant to this pull request? This may include screenshots, references, and/or implementation notes.

No.

Checklist

Please ensure the following tasks are completed before submitting this pull request.

AI Assistance

When authoring the changes proposed in this PR, did you use any kind of AI assistance?

  • Yes
  • No

If you answered "yes" above, how did you use AI assistance?

  • Code generation (e.g., when writing an implementation or fixing a bug)
  • Test/benchmark generation
  • Documentation (including examples)
  • Research and understanding

Disclosure

If you answered "yes" to using AI assistance, please provide a short disclosure indicating how you used AI assistance. This helps reviewers determine how much scrutiny to apply when reviewing your contribution. Example disclosures: "This PR was written primarily by Claude Code." or "I consulted ChatGPT to understand the codebase, but the proposed changes were fully authored manually by myself.".

This PR was written primarily by Claude Code. I found the bug by reading the stats/base/dists quantile functions for 1 - p cancellation.


@stdlib-js/reviewers

`1 - p` rounds away the low-order bits of `p`. For `p < 2^-53` the
function returned `-1`, outside the support, and for small `lambda` it
was off by more than 10%.
@Abhist17
Abhist17 requested a review from a team October 1, 2026 23:19
@stdlib-bot stdlib-bot added Statistics Issue or pull request related to statistical functionality. Needs Review A pull request which needs code review. labels Oct 1, 2026
@stdlib-bot

stdlib-bot commented Oct 1, 2026 •

Copy link
Copy Markdown
Contributor

Coverage Report

Package Statements Branches Functions Lines
stats/base/dists/planck/quantile $\\color{green}288/288$
$\\color{green}+0.00\\%$
$\\color{green}25/25$
$\\color{green}+0.00\\%$
$\\color{green}4/4$
$\\color{green}+0.00\\%$
$\\color{green}288/288$
$\\color{green}+0.00\\%$

The above coverage report was generated for the changes in this PR.

return PINF;
}
return ceil( -ln( 1.0-p ) / lambda ) - 1.0;
return ceil( -log1p( -p ) / lambda ) - 1.0;

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Suggested change
return ceil( -log1p( -p ) / lambda ) - 1.0;
return max( 0.0, ceil( -log1p( -p ) / lambda ) - 1.0 );

Should this clamp at 0 like geometric/quantile? -log1p( -p ) / lambda can still underflow to 0, so quantile( 5.0e-324, 2.0 ) and quantile( 1.0e-300, 1.0e300 ) return -1. The clamp would also settle the p = 0 question in the description. Same in lib/factory.js:73 and src/main.c:43 (plus @stdlib/math/base/special/max in manifest.json).

Copy link
Copy Markdown
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Yes, clamped with max( 0.0, ... ) in main.js, factory.js and src/main.c, with @stdlib/math/base/special/max added to all three manifest confs (571b8e9). New tests in all three files cover p = 0, quantile( 5.0e-324, 2.0 ) and quantile( 1.0e-300, 1.0e300 ), which returned -1 before and now return 0. That also answers the p = 0 question, so I dropped it from the description.

t.end();
});

tape( 'the created function evaluates the quantile for small `p`', function test( t ) {

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Suggested change
tape( 'the created function evaluates the quantile for small `p`', function test( t ) {
tape( 'the returned function evaluates the quantile for small `p`', function test( t ) {

Match the other test names in this file.

Copy link
Copy Markdown
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Renamed in 571b8e9.

`-log1p( -p ) / lambda` can still underflow to `0` (e.g. `p = 5e-324`,
`lambda = 2`), giving `ceil( 0 ) - 1 = -1`, outside the support. Clamping
at `0` like `geometric/quantile` fixes that and also returns `0` at `p = 0`.

Signed-off-by: Abhist17 <abhistkamle.170903@gmail.com>

This branch has not been deployed

No deployments
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

Needs Review A pull request which needs code review. Statistics Issue or pull request related to statistical functionality.

Projects

None yet

Development

Successfully merging this pull request may close these issues.

3 participants