Repository navigation
Conversation
…alars With fvar<var> or fvar<fvar<var>> coefficients the GLM lpmfs compute their partials in var / fvar<var>, and the Eigen products x * beta and x^T * theta_derivative then run as coefficient loops that put two autodiff nodes per matrix entry on the tape (2 n K each). Calling multiply instead resolves to the reverse-mode overload (one node, dense product inside the callback) or the forward-mode overload, so one second-order sweep of a GLM log density adds O(n) nodes instead of O(n K). Measured on an embedded-Laplace hierarchical GLM (n = 30..400 rows, 49 coefficients): 214 -> 22 nodes per row per second-order sweep, 508 -> 131 per third-order sweep, laplace_marginal 2.5-3.3x faster; values and gradients unchanged at every autodiff order. categorical_logit_glm_lpmf with a matrix beta did not compile with forward-mode scalars before (mixed-scalar matrix-matrix product in Eigen's GEMM kernel); the products on the x-partials path go through multiply as well, and the mix test gains a three-category case. Co-Authored-By: Claude Fable 5.1 <[email protected]>
…logit_glm_lpmf multiply of a compile-time row vector and a compile-time column vector is a dot product returning a scalar, so the single-row broadcast case combined with a one-column beta (a vector) had nothing to call rowwise() on. That branch is cheap and keeps the Eigen product, as the single-row branches of the other GLM densities do. Co-Authored-By: Claude Fable 5.1 <[email protected]>
|
Not sure why the CI is failing, I think jenkins is running out of memory or something like that. |
|
Yeah we just moved our jenkins over to a new set of servers recently so there are some hiccups I'm looking at it now |
|
I think the meta here is that we really want |
|
Yes, I agree, this feels a bit like hacking and it's quite non-obvious to a new developer why this is even a semantic difference in the first place. If there was a more conceptual way to improve this that would be great, but it's probably outside my ability to address at this point. I just noticed this is a way I could speed-up an application I'm building on top of stan that has a hierarchical GLM model, so this has a practical motivation for me and I went for the approach I could do myself. |
|
Yes tbc this is a good PR and once I restart the jenkins I think this is good to merge as is. We have not done the matrix multiply overloads because it requires messing with a lot of Eigens internals. I'm taking a crack at that, but for now I think what you have here is best |
|
I think this is not jenkins being flaky, but legitimately running out of memory This PR touches 9 GLM headers, so I think this needs a restart with higher memory available. Is either of that possible? |
|
I've temporarily increased the memory available to that stage to test your theory |
|
Thanks @WardBrian, but seems like it wasn't enough. I triggered the rerun by merging develop. I still crashes while compiling, the tests don't even get to run. I could perhaps try breaking this up into up to 7 PRs, but that feels like spamming. Adding even more memory or setting Thoughts? |
|
Sorry, the way I did it was to trigger https://jenkins.flatironinstitute.org/job/CCM/job/Stan/job/math/job/PR-3418/9/ with the temporary increase, which then got killed by your push. https://jenkins.flatironinstitute.org/job/CCM/job/Stan/job/math/job/PR-3418/11/ is running now |
|
@jachymb it looks like the quick tests now passed, and it is the expression tests that are failing (I believe legitimately?) |
With a column-vector x (one attribute), multiply(x^T, theta_derivative) is a dot product and its scalar result cannot be assigned to the beta partials; this broke the expression tests of normal_id_glm_lpdf, whose Stan signatures take x as a vector. Store dot_product(...) in the single partial. In categorical_logit_glm_lpmf a row-vector x gives the same shape with a column-vector beta (linear predictor) and a row-vector beta (x partial).
|
Hi, @WardBrian, thanks for the help. Yes it was a legitimate fail, apparently I didn't run the expression tests first on my machine. Sorry about that. Anyway, it should be fixed now. The issue was with May I ask for a re-trigger with the increased memory? I re-triggered it with the fix push, but got the OOM again. |
Jenkins Console Log Machine informationDistributor ID: Ubuntu Description: Ubuntu 20.04.3 LTS Release: 20.04 Codename: focal CPU: Architecture: x86_64 CPU op-mode(s): 32-bit, 64-bit Byte Order: Little Endian Address sizes: 52 bits physical, 57 bits virtual CPU(s): 192 On-line CPU(s) list: 0-191 Thread(s) per core: 2 Core(s) per socket: 48 Socket(s): 2 NUMA node(s): 2 Vendor ID: AuthenticAMD CPU family: 25 Model: 17 Model name: AMD EPYC 9474F 48-Core Processor Stepping: 1 Frequency boost: enabled CPU MHz: 1494.432 CPU max MHz: 4114.4229 CPU min MHz: 1500.0000 BogoMIPS: 7189.39 Virtualization: AMD-V L1d cache: 3 MiB L1i cache: 3 MiB L2 cache: 96 MiB L3 cache: 512 MiB NUMA node0 CPU(s): 0-47,96-143 NUMA node1 CPU(s): 48-95,144-191 Vulnerability Gather data sampling: Not affected Vulnerability Indirect target selection: Not affected Vulnerability Itlb multihit: Not affected Vulnerability L1tf: Not affected Vulnerability Mds: Not affected Vulnerability Meltdown: Not affected Vulnerability Mmio stale data: Not affected Vulnerability Old microcode: Not affected Vulnerability Reg file data sampling: Not affected Vulnerability Retbleed: Not affected Vulnerability Spec rstack overflow: Mitigation; Safe RET Vulnerability Spec store bypass: Mitigation; Speculative Store Bypass disabled via prctl Vulnerability Spectre v1: Mitigation; usercopy/swapgs barriers and __user pointer sanitization Vulnerability Spectre v2: Mitigation; Enhanced / Automatic IBRS; IBPB conditional; STIBP always-on; PBRSB-eIBRS Not affected; BHI Not affected Vulnerability Srbds: Not affected Vulnerability Tsa: Mitigation; Clear CPU buffers Vulnerability Tsx async abort: Not affected Vulnerability Vmscape: Mitigation; IBPB before exit to userspace Flags: fpu vme de pse tsc msr pae mce cx8 apic sep mtrr pge mca cmov pat pse36 clflush mmx fxsr sse sse2 ht syscall nx mmxext fxsr_opt pdpe1gb rdtscp lm constant_tsc rep_good amd_lbr_v2 nopl xtopology nonstop_tsc cpuid extd_apicid aperfmperf rapl pni pclmulqdq monitor ssse3 fma cx16 pcid sse4_1 sse4_2 x2apic movbe popcnt aes xsave avx f16c rdrand lahf_lm cmp_legacy svm extapic cr8_legacy abm sse4a misalignsse 3dnowprefetch osvw ibs skinit wdt tce topoext perfctr_core perfctr_nb bpext perfctr_llc mwaitx cpb cat_l3 cdp_l3 hw_pstate ssbd mba perfmon_v2 ibrs ibpb stibp ibrs_enhanced vmmcall fsgsbase bmi1 avx2 smep bmi2 erms invpcid cqm rdt_a avx512f avx512dq rdseed adx smap avx512ifma clflushopt clwb avx512cd sha_ni avx512bw avx512vl xsaveopt xsavec xgetbv1 xsaves cqm_llc cqm_occup_llc cqm_mbm_total cqm_mbm_local user_shstk avx512_bf16 clzero irperf xsaveerptr rdpru wbnoinvd amd_ppin cppc arat npt lbrv svm_lock nrip_save tsc_scale vmcb_clean flushbyasid decodeassists pausefilter pfthreshold avic v_vmsave_vmload vgif x2avic v_spec_ctrl vnmi avx512vbmi umip pku ospke avx512_vbmi2 gfni vaes vpclmulqdq avx512_vnni avx512_bitalg avx512_vpopcntdq la57 rdpid overflow_recov succor smca fsrm flush_l1d debug_swap G++: g++ (Ubuntu 9.4.0-1ubuntu1~20.04) 9.4.0 Copyright (C) 2019 Free Software Foundation, Inc. This is free software; see the source for copying conditions. There is NO warranty; not even for MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. Clang: clang version 10.0.0-4ubuntu1 Target: x86_64-pc-linux-gnu Thread model: posix InstalledDir: /usr/bin |
|
Seems like it passed now. I recommend this gets merged before other things so we don't have to rerun it and deal with the OOM again after any |
This branch has not been deployed
4 participants
AI use disclosure: The code as well as this commentary was done with the help of claude Fable 5.1
Summary
The GLM log densities compute their two design-matrix products,
x * betaandx^T * theta_derivative, withmultiplyinstead of Eigen'soperator*. With forward-mode coefficients (Hessian-vector products, the third-order term of the embedded Laplace approximation) this puts one autodiff node per result coefficient on the tape instead of two per matrix entry, which makes such sweeps several times faster.Why it matters: the partials are computed in
partials_return_t, which isvarforfvar<var>inputs, and Eigen evaluates adouble-by-varproduct coefficient by coefficient,2 n Knodes per product on every sweep.multiplyis resolved by argument-dependent lookup at instantiation (reverse-mode overload forvarpartials, forward-mode overload forfvar<var>partials), soprimgains no dependency onrev. Fordoubleandvarinputs theprimoverload returns the same lazy product as before.categorical_logit_glm_lpmfwith a matrixbetadid not compile withfvar<var>at all (Eigen's GEMM kernel cannot mix scalar types), which the existing one-column test never exercised; itsx-partials products are switched as well.Measured on an embedded-Laplace hierarchical GLM (49 coefficients, 30 to 400 observations per group): 214 to 22 nodes per observation per second-order sweep, 508 to 131 per third-order sweep,
laplace_marginal3.1 to 3.8 times faster. Values and gradients match develop to 12 digits at every autodiff order.Tests
test/unit/math/mix/prob/glm_forward_mode_tape_test.cpp(new): onefvar<var>sweep of five GLM densities must add fewer than40 nnodes to the tape; before the change it added about4 n K + 20 n.categorical_logit_glm_lpmf_test.cpp: new three-categoryexpect_adcase, which did not compile before.primandmixGLM tests are unchanged.Side Effects
None for
doubleandvararguments. With forward-mode arguments the reverse-modemultiplycopies the design matrix into the arena once per call (n Kdoubles), which is what replaces the2 n Knode allocations.Release notes
GLM log densities evaluated with forward-mode scalars (nested autodiff, for example the embedded Laplace approximation) now create O(n) instead of O(n K) autodiff nodes per sweep.
categorical_logit_glm_lpmfnow compiles with forward-mode coefficients for more than one category.Checklist
Copyright holder: Jachym Barvinek
The copyright holder is typically you or your assignee, such as a university or company. By submitting this pull request, the copyright holder is agreeing to the license the submitted work under the following licenses:
- Code: BSD 3-clause (https://opensource.org/licenses/BSD-3-Clause)
- Documentation: CC-BY 4.0 (https://creativecommons.org/licenses/by/4.0/)
the basic tests are passing
./runTests.py test/unit)make test-headers)make test-math-dependencies)make doxygen)make cpplint)the code is written in idiomatic C++ and changes are documented in the doxygen
the new changes are tested