Performance follow-up: GGM chol deferral, NUTS gradient workspaces, mixed AM caching, bgmCompare trims#166
Merged
Conversation
Codecov Report❌ Patch coverage is
Additional details and impacted files@@ Coverage Diff @@
## develop #166 +/- ##
===========================================
+ Coverage 90.95% 91.00% +0.05%
===========================================
Files 104 104
Lines 14900 14914 +14
===========================================
+ Hits 13552 13573 +21
+ Misses 1348 1341 -7 ☔ View full report in Codecov by Harness. 🚀 New features to boost your workflow:
|
…ement The Metropolis, tuning, and edge-indicator sweeps read only Sigma and log|K| between accepts, so the per-accept Givens update/downdate pair (2 x O(p^2)) and the p-log get_log_det refresh are dropped. Accepted edge/diagonal moves advance log_det_precision_ by the same determinant-lemma ratio the MH ratio uses (O(1) from cached Sigma), and every sweep ends with one refresh_cholesky(), which rebuilds the factor and resets Sigma and the log-det exactly. This replaces the per-sweep drift-tolerance check: the SMW floating-point drift is now cleared unconditionally at sweep end. Verification: test-ggm-gibbs slow gates 20/20, test-hier-zratio-identity 14/14 (zero chol symmetry warnings), test-sample-ggm-prior 53/53, test-ggm-nuts fast 2/2.
The Z-ratio tilt-grid characteristic function and the SBM cluster draw were the last per-sweep arma::exp sites outside the platform macros. No change on macOS/Linux; Windows builds avoid UCRT exp here. Verification: zratio-engine and sbm-correction suites green.
…ain component log_pseudoposterior_main_component bound the per-group count and Blume-Capel sufficient-statistic matrices by value on every evaluation; bind them by const reference, matching the other component functions. Replace the zero-filled V x max_categories and V x V scratch matrices, of which only one row and one column were used, with the group main-effect vector already computed and a length-V pairwise-weight vector hoisted out of the group loop. Verified: tests/testthat/test-bgmCompare.R passes (40 passed, 0 failed).
…-data imputation impute_missing_bgmcompare rebuilt the V x V group pairwise-effect matrix per missing cell although the parameters are constant across the call; build it once per group before the imputation loop. Mutate the per-group count and Blume-Capel statistics in place instead of copying the matrix out and back per changed cell. Replace the full O(n_g * V^2) crossproduct rebuild of the group pairwise statistics with the O(V) rank-1 update of row and column `variable`, with the diagonal corrected by delta^2. Verified: tests/testthat/test-bgmCompare.R passes with NOT_CRAN=true (42 passed, 0 failed); a seeded bgmCompare fit with na_action = "impute" (adaptive-metropolis, 2 chains) produces raw_samples identical() to the pre-change build.
…pair component log_pseudoposterior_pair_component converted two integer observation columns to double on every evaluation. Both adaptive-Metropolis pairwise sweeps already hold per-group double observation matrices (obs_double_groups); pass those down and index the needed column directly instead of converting per call. Verified: tests/testthat/test-bgmCompare.R passes with NOT_CRAN=true (42 passed, 0 failed).
log_ggm_ratio_edge/_diag evaluate the quadratic-form difference through the low-rank structure of the proposed covariance change (O(nq + q^2) dot products) instead of rebuilding the proposed conditional mean via an O(nq^2) GEMM and both n x q quadratic forms from scratch. The mean and cross moves (update_continuous_mean, update_pairwise_cross, update_edge_indicator_cross) evaluate the GGM term as the cached value plus a rank-1 quadratic-form delta, so conditional_mean_ is only touched on accept and the per-proposal n x q save/restore copies are gone. cond_mean_scratch_ is removed with its last consumers. Verification: fast mixed suites pass (test-mixed-mrf-likelihood 30, test-mixed-gradient 8, test-mixed-edge-indicators 3, test-mixed-theta-space 18). Same-seed 8+8 AM regression (iter 52, with/without selection): indicator paths bitwise identical; draws agree to 7.8e-16 over a short horizon and to 1.8e-7 over 104 iterations (Robbins-Monro feedback amplifies last-bit ln_alpha rounding differences; no acceptance decision flips). AM 8+8 selection fit (iter 1000, load_all -O0): 20.5 s -> 17.5 s.
Accepted Kyy (off-diagonal, diagonal, edge-indicator) and cross (parameter, edge-indicator) moves adopt the caches the ratio evaluation already computed instead of rebuilding them via recompute_marginal_interactions + recompute_am_caches (an O(np^2) GEMM, all p OMRF denominators, and log_conditional_ggm per accept): - adopt_kyy_proposal_caches swaps marginal_matvec_/ll_marginal_cache_ with the proposal scratch, applies 2 A_xy dSigma A_xy' to marginal_interactions_/cross_term_ via the E = dSigma A_xy' buffer cached by omrf_ratio_for_covariance_change, shifts conditional_mean_ by the low-rank factors of the accepted ratio, and advances ll_ggm_cache_ by the ratio value. - adopt_cross_proposal_caches applies the rank-2 row/column-i update to marginal_interactions_/cross_term_, rank-1-updates cross_matvec_.col(j) and conditional_mean_, and copies the proposed bias and marginals. - The diagonal of marginal_interactions_ is overwritten with mdiag_prop_ so it matches the adopted per-variable marginals. - cholesky_update_after_precision_edge/_diag report whether the factors advanced incrementally; the rebuild fallback keeps the full recompute_am_caches, since the proposal buffers no longer match a from-scratch decomposition. - recompute_am_caches now refreshes marginal_interactions_ first, so the sweep-top full refresh (do_one_metropolis_step, sweep_within_model_mh, update_edge_indicators) bounds the floating-point drift of the low-rank cache updates to one sweep. Verification: fast mixed suites pass (30/8/3/18 assertions). Same-seed 8+8 AM regression (iter 52, with/without selection): indicator paths bitwise identical; draws agree to 6.5e-16 over a short horizon and to 2.4e-7 over 104 iterations (Robbins-Monro amplification of last-bit ln_alpha rounding; no acceptance decision flips). AM 8+8 selection fit (iter 1000, load_all -O0): 17.5 s -> 13.5 s (20.5 s at the audit base).
…rd pass forward_map allocated its entire result per call — Phi, K, and per constrained column the Q/R QR matrices plus the A_q build and .t() temporaries — and theta_gradient_from_phi_bar allocated four q x q / q x m_q matrices per column, totalling O(p^3) doubles per leapfrog step on sparse graphs. The engine now owns the ForwardMapResult and the backward-pass scratch as mutable members reused across calls (zeros(), eye(), and copy assignment keep the existing allocation when dimensions match): forward_map returns a const reference valid until the next call, A_q^T is written directly into the in-place QR working matrix, and P / Phi_bar live in the workspace. Each chain clone owns its engine copy and the sampler is single-threaded within a chain. Verification: same-seed GGM NUTS + selection and mixed NUTS + selection draws are bitwise identical to the previous build (identical() on the full raw samples); ggm-gradient, ggm-nuts, mixed-gradient, mixed-theta-space, and ggm-gibbs suites green.
The OMRF gradient issued one X^T * E GEMV per variable and the mixed gradient two (pw_grad and diff_pw), plus an O(p*q) scalar loop for the cross-contribution. Collect the per-variable expected-score vectors into an n x p buffer and assemble the pairwise, cross, and Theta-bar contributions after the loop from batched GEMMs: one X^T * E_all in the OMRF, two in the mixed gradient, with the cross self- and cross-contributions combined as 4 (diff_pw_all + diff_pw_all^T) * (A_xy Sigma_yy) and the mean gradient as one GEMV. Verification: fixed-seed ordinal NUTS draws (p=10, n=300) are bitwise identical; mixed 6+6 NUTS draws deviate by max 3.1e-10 over a 10+10-iteration horizon (reassociated BLAS accumulation). FD gradient suites pass: test-mixed-gradient.R (8), test-mixed-theta-space.R (18), test-ggm-gradient.R (14).
The numerator loop recomputed the per-variable rest score verbatim from X_marginal — a second full O(np) pass with p n-vector allocations per leapfrog step. Accumulate the quadratic self-term, dot(x_s, rest), and the main-effect sufficient-statistic sums inside the Part-1 loop, which already holds the rest vector, and drop the second pass. Verification: fixed-seed ordinal NUTS draws remain bitwise identical; mixed 6+6 NUTS draws deviate by max 1.3e-9 from the pre-branch baseline over a 10+10-iteration horizon (logp accumulation order). FD gradient suites pass: test-mixed-gradient.R (8), test-mixed-theta-space.R (18), test-ggm-gradient.R (14).
MaartenMarsman
force-pushed
the
perf/audit-2026-07-09-follow-up
branch
from
6e86700 to
51701f7
Compare
MaartenMarsman
changed the title
Collaborator
Author
|
Revised per review; the branch is rebuilt on develop @ dac1627 (force-pushed, now 10 commits). Dropped: chain_backend fork dispatcher. Removed entirely — R interface, validation, tests, docs. Dropped: prior-density inlining (was 2c7d94e). An isolated A/B first showed -5.7% on the ordinal NUTS cell, but a control rebuild of functionally identical code moved that same cell by 3.3% — build-layout noise. With that floor established, neither the ordinal nor the GGM cell (which loses ~7M Rmath calls) shows a gain that survives. No measurable benefit, so the maintenance cost isn't justified. Dropped: NUTS leapfrog memoizer (was 85091d4). Instrumented counters confirmed the mechanism is real but depth-dependent: it saves one joint evaluation per direction flip, i.e. a fraction of roughly (depth-1)/2^(depth+1) of all evaluations. That's 7.8% on the shallow GGM cell (matching the measured 7.7% wall-clock) but ~2% once warmup trees get deep, as on the ordinal cell. Dropped per review; the ordinal NUTS row in the table loses its earlier apparent gain, which the noise-floor control attributes to layout noise anyway. Fixed: unused pairwise_g in the imputation commit (now 3a8c8eb). That matrix was a leftover alias — the code built per-group effect matrices but also kept a redundant copy. The commit now writes the per-group matrices directly and pairwise_g is gone; imputation draws remain bitwise identical to develop. Re-verified on the rebuilt branch: full suite with NOT_CRAN 7000/0/0, R CMD check --as-cran 0 errors / 0 warnings / 2 notes, benchmarks rerun from fresh -O2 installs (table in the description updated). |
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.
Implements the remaining items from the 2026-07-09 performance re-audit (dev/audit/2026-07-05-perf-fixes-status.md), 10 commits grouped by area.
GGM (AM/Gibbs)
NUTS
Mixed AM
bgmCompare
Platform
Verification
Benchmarks (installed -O2 builds vs develop @ dac1627)
The ordinal NUTS cell is unchanged: functionally identical rebuilds of this package differ by up to ~3% on that cell (build-layout noise), and the residual difference is inside that floor. Its gradient is already dominated by the category loop, which is untouched.