Repository navigation
Speed up unstratified AE-term confidence intervals (~18x) - #257
Conversation
biroot() called the objective func_d() once per grid point -- ~100 times per confidence interval, and once per AE term -- so the scan was ~1M interpreted scalar calls on a large forestly example, dominating extend_ae_specific_inference(). func_d() is already written with vectorized arithmetic. In the unstratified case its argument `d` can be a vector, so evaluate the whole scan grid in one call, then bisect only the intervals whose endpoints change sign. The stratified case (which reduces over strata with sum() and cannot take a vector `d`) still evaluates the grid point by point. func_d()'s unstratified special case is rewritten as a vectorized mask. Results are unchanged up to floating-point tolerance: max abs difference 4.4e-16 across 4,000 unstratified and 500 stratified random cases. On a 9,510-term forestly example prepare drops from ~132s to ~39s. Co-Authored-By: Claude Opus 4.8 <noreply@anthropic.com>
Even with the vectorized grid scan, extend_ae_specific_inference() called rate_compare_sum() once per AE term -- thousands of scans, each building the scalar machinery and running its own bisection. That per-term call overhead dominated prepare on large analyses. Add rate_compare_sum_batch(), a vectorized-across-terms reimplementation of rate_compare_sum()'s unstratified path: it evaluates the bisection grid for all terms at once (an nt x E matrix) and refines every sign-change bracket with a single vectorized bisection loop. When there is no stratification, extend_ae_specific_inference() now computes all terms' confidence intervals with one call instead of a per-term loop; stratified inputs still use rate_compare_sum(). Value-preserving: max abs difference 0 (bit-identical est/z_score/p/ lower/upper) against rate_compare_sum() across all 9,510 terms of a large forestly example, plus random and edge-case (zero-event, tiny n, nonzero delta, two-sided) checks. On that example prepare_ae_forestly() drops from ~39s to ~2s, on top of the earlier grid vectorization. Co-Authored-By: Claude Opus 4.8 <noreply@anthropic.com>
yihui
changed the title
The constrained-MLE variance algebra (the l0..l3 -> q -> p -> r0t -> r1t -> vart block) was copy-pasted four times: the score statistic and CI objective in rate_compare_sum(), and again in rate_compare_sum_batch()'s func_pts() and point-estimate block. Pull it into a single vectorized mn_vart() helper that every caller shares; each op is elementwise, so the CI grid scan (scalar term params, vector d) and the across-terms batch (equal-length vectors) both use it via recycling. An adjust_p flag reproduces each caller's exact prior behavior: the score statistic path never nudged p off zero, the CI objective did. The batch point-estimate previously applied the nudge but its z_score matched the scalar path bit-for-bit on all test data (the nudge never fired); it now uses adjust_p = FALSE to match the scalar path by construction. Results unchanged: scalar (one/two-sided, delta != 0), stratified, and batch-vs-scalar all match the pre-refactor code to maxdiff 0. Co-Authored-By: Claude Opus 4.8 <noreply@anthropic.com>
If it is stratified, then will we have the same running speed as before? |
…ce tests Review follow-up on #257: - Rename rate_compare_sum_batch() to rate_compare_sum_unstratified() (and the local use_batch -> use_unstratified, roxygen title) so the name reads in domain terms rather than algorithm-engineering jargon. - Move the chisq_crit and unstratified definitions above biroot() so `unstratified` is defined before biroot() reads it, instead of appearing below the function that uses it. - Add test-independent-testing-rate_compare_sum_unstratified.R: assert the vectorized unstratified path matches looping the scalar rate_compare_sum() one term at a time, to 1e-10 (est/z_score/p) and 1e-8 (CI limits), across random inputs, edge cases (zero/all events, single subject, extreme split), nonzero delta + two-sided, and NA handling. rate_compare() and biroot() are exercised through rate_compare_sum(); extend_ae_specific()'s unstratified path runs through this function and is covered by the existing extend_ae_specific tests. Co-Authored-By: Claude Opus 4.8 <noreply@anthropic.com>
Essentially yes — stratified inputs run at about the same speed as before this PR. The two big wins here (evaluating the whole bisection grid in one vectorized call, and computing all AE terms together in They do keep the path-agnostic improvements (hoisted loop invariants, each grid edge evaluated once and shared between adjacent intervals, and the shared |
Review follow-up on #257: - Add tests/testthat/helper-rate_compare_sum_old.R: a verbatim frozen copy of rate_compare_sum() from 26ba273^ (before this PR vectorized the grid scan and added the across-terms path), renamed rate_compare_sum_old(). The equivalence tests now compare BOTH the current scalar rate_compare_sum() and the new rate_compare_sum_unstratified() against this genuine original, instead of against each other. - Run every case at alpha = 0.025 and 0.05. - Drop the explicit weight = "ss" from the reference loop. For a single unstratified term the stratum weight normalizes to 1, so est/z_score/p and the CI do not depend on weight; a comment notes this. Co-Authored-By: Claude Opus 4.8 <noreply@anthropic.com>
LittleBeannie
left a comment
LittleBeannie
left a comment
There was a problem hiding this comment.
This is an amazing PR! I can't imagine how much thought, intelligence, and effort went into this.
|
The hard work was mostly done by AI. The PRs on the forestly side are more demanding in effort, but I think I've finished the hardest bones now. I just need to review all |
2 participants
Summary
Follow-up to #253. Computing the Miettinen-Nurminen risk-difference confidence intervals in
extend_ae_specific_inference()was the dominant cost of preparing an interactive AE forest plot on large analyses. Two layers of overhead remained:biroot()called the objectivefunc_d()once per grid point (~100× per CI).extend_ae_specific_inference()calledrate_compare_sum()once per AE term — thousands of independent scans, each rebuilding the scalar machinery and running its own bisection.This PR removes both, for the unstratified case (the forest-plot path):
biroot()evaluates the whole scan grid in one vectorizedfunc_d()call, then bisects only the intervals whose endpoints change sign. Stratified inputs (which reduce over strata withsum()and can't take a vectord) still evaluate the grid point by point.rate_compare_sum_batch()reimplementsrate_compare_sum()'s unstratified path over vectors of(n0, n1, x0, x1). It evaluates the grid for all terms at once (annt × Ematrix) and refines every sign-change bracket with a single vectorized bisection loop.extend_ae_specific_inference()uses it whenever there is no stratification; stratified inputs still userate_compare_sum()per term. The publicrate_compare_sum()API is unchanged.Correctness
Value-preserving. Against
rate_compare_sum():0) onest,z_score,p,lower,upper.(n0, n1, x0, x1)(one-sided, two-sided,delta = 0.1) and edge cases (zero events, tinyn,delta = ±): max abs diff0.End-to-end, the prepared
ci_lower/ci_uppermatch a direct scalarrate_compare_sum()computation exactly.(One behavior change: the batch path does not emit the per-term "no CI limit found"
message()that the scalar path prints for terms with no root in range; the returnedNAlimits are the same.)Impact
On a 9,510-term forestly example,
prepare_ae_forestly():~60× overall; ~18× beyond #253.
🤖 Generated with Claude Code