Add fastpidc: standalone Python port of FastPIDC.jl - #10
Merged
Merged
Conversation
Adds a NumPy/SciPy implementation of the MI, CLR, PUC and PIDC network inference algorithms under python/, managed with uv, that does not call out to Julia at runtime. Ported module-for-module (discretizers, probability estimators/formulae, network inference, I/O) and validated against a freshly run FastPIDC.jl (see python/tests/test_against_julia.py) to float64 precision on both the CPU and GPU backends. GPU acceleration is implemented in plain CUDA C (python/src/fastpidc/kernels/pidc_kernels.cu), ported from and numerically verified against FastPIDCCUDAExt's chunked kernels, so the two packages share one canonical kernel source; Python's cuda backend loads it via cupy today, with wiring the Julia extension to the same file left as documented follow-up (see the file's header comment). While cross-validating, found a narrow off-by-one bug in FastPIDC.jl's Bayesian-blocks change-point backtracking (a BoundsError, silently caught and downgraded to a uniform-width fallback) for data with very few unique values whose optimal segmentation assigns every point its own block. The Python port's 0-indexed rewrite does not reproduce it; see the comment in python/tests/test_against_julia.py::test_all_algorithms_match_julia. Also documents the two-package repository layout in the root README, and adds a path-scoped Python CI workflow alongside the existing Julia one. Co-Authored-By: Claude Sonnet 5 <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_01EAVcaxbtkiPMc4PQBmb4Tk
…-blocks off-by-one FastPIDCCUDAExt now compiles python/src/fastpidc/kernels/pidc_kernels.cu with nvcc (targeting the active device's compute capability) and drives it via CUDA.jl's CuModule/cudacall, replacing its own CUDA.jl-native kernel copy. Device buffers use Julia's column-major layout with dimensions reversed relative to the kernel's documented row-major shapes, so the flat in-memory layout matches exactly and no transposition is needed at the call boundary; bin ids (1-indexed in Julia) are shifted down by one before upload. Verified against the CPU backend (float64-precision match on yeast1_10 and toy_small_200) and against the Python cuda backend running the identical .cu file (2e-15 max diff), and the full Julia test suite passes, including the 1000-gene CUDA benchmark (208x speedup vs. 10-worker CPU, matching to 3e-12). Also fixes the off-by-one bug identified in PR #10: binedges_bayesian_blocks' change-point backtracking under-allocated its `change_points` array by one slot, causing a BoundsError (caught and silently downgraded to a uniform-width fallback) whenever the optimal segmentation assigned every point its own block - e.g. mostly-constant data with a single outlier. Sized the array for the true worst case instead. With this fixed, the Python port's bin counts now match FastPIDC.jl's exactly on every gene in the test fixtures (previously a handful legitimately diverged); the cross-validation test is updated to assert that directly instead of carving out exceptions. It also surfaced a separate, non-bug numerical sensitivity: CLR/PIDC's (score - mean) / sqrt(variance) context weighting is inherently ill-conditioned for a gene whose MI against everything is zero up to floating-point noise (background variance ~1e-30), so ordinary summation-order differences between the two implementations can move the standardized score by O(1). MI and PUC (unnormalized) still compare bit-for-bit; CLR/PIDC now compare by rank correlation instead of exact equality, with the reasoning documented in the test. Also adds `LinearAlgebra` to Project.toml's test targets (needed by test/cuda_numeric_tests.jl's `issymmetric`, missing before this change and only surfaced when running `Pkg.test()` in its sandboxed environment) and regenerates the stale, untested test/baseline_outputs/pidc_toy_edges.tsv snapshot now that the discretizer fix changes its values. Co-Authored-By: Claude Sonnet 5 <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_01EAVcaxbtkiPMc4PQBmb4Tk
Member
Author
|
@jjschirle I'm not sure what it means to "review" an entirely new package, but here it is. Let me know if you see any potential issues. |
Resolves a conflict in the CUDA extension between two independent changes to compute_puc_full_cuda: wiring device buffers to the shared kernel source's reversed-dims layout (this branch) and sizing chunks to the currently-free GPU memory instead of a fixed 256 (master). Keeps both: reversed-dims buffer shapes with memory-adaptive chunk sizing. Co-Authored-By: Claude Sonnet 5 <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_01EAVcaxbtkiPMc4PQBmb4Tk
python-CI.yml never instantiated the root Julia project, so on ubuntu-latest runners (which ship a preinstalled julia) the Julia cross-validation tests failed instead of skipping. Also fixes docs/make.jl failures: the CUDA extension's @docs block referenced a nonexistent binding (bayesian_blocks_dp_kernel! is a CUDA C kernel name, not a Julia function), and solve_bayesian_blocks_cuda/bayesian_blocks_cuda_available were missing docstrings despite being @ref'd. Co-Authored-By: Claude Sonnet 5 <noreply@anthropic.com>
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
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.
Summary
python/: fastpidc, a standalone NumPy/SciPy port of the MI, CLR, PUC and PIDC network inference algorithms, managed with uv. It does not call out to Julia at runtime.FastPIDC.jl(python/tests/test_against_julia.py, skipped automatically whenjuliaisn't onPATH) to float64 precision, on both CPU and GPU backends.python/src/fastpidc/kernels/pidc_kernels.cu). Python'scudabackend compiles it viacupy/nvrtc;FastPIDCCUDAExt(Julia) now compiles the same file withnvccand drives it withCUDA.jl'sCuModule/cudacall, replacing its previous CUDA.jl-native kernel copy. Verified end-to-end on a real GPU: Julia-CUDA vs Julia-CPU, Python-CUDA vs Python-CPU, and Julia-CUDA vs Python-CUDA all agree to float64 precision, and the full Julia test suite passes (including the 1000-gene benchmark: 208x speedup vs. 10-worker CPU).README.md, and adds a path-scoped.github/workflows/python-CI.ymlalongside the existing JuliaCI.yml, so each package's CI is independent.Bug fixes in the Julia source
src/discretizers.jl): the change-point backtracking under-allocated itschange_pointsarray by one slot, raising aBoundsError(silently caught and downgraded to auniform_widthfallback) whenever the optimal segmentation assigned every point its own block - e.g. mostly-constant data with a single outlier (repro:FastPIDC.binedges(FastPIDC.DiscretizeBayesianBlocks(), [fill(0.0, 999); 0.5])). Fixed by sizing the array for the true worst case. With this fixed, the Python port's bin counts now match FastPIDC.jl's exactly on every gene in the test fixtures.(score - mean) / sqrt(variance)context weighting is inherently ill-conditioned for a gene whose MI against everything is zero up to floating-point noise (background variance ~1e-30), so ordinary summation-order differences between the two implementations can move the standardized score by O(1). MI and PUC (unnormalized) compare bit-for-bit between Julia and Python; CLR/PIDC now compare by rank correlation instead of exact equality in the cross-validation test, with the reasoning documented there.test/cuda_numeric_tests.jlwas missingusing LinearAlgebraforissymmetric, which only surfaces when runningPkg.test()in its sandboxed test environment (addedLinearAlgebratoProject.toml's test targets too).Test plan
cd python && uv run pytest -q— 69 passed (includes live cross-checks against a freshly runjuliaprocess and against a real GPU viacupy, both available in this sandbox)cd python && uv run ruff check ./ruff format --check .— cleanjulia --project=. -e 'import Pkg; Pkg.test()'— full Julia suite passes on a real GPU, including CUDA numeric equivalence and the 1000-gene benchmarkpython-CI.yml) will run on a GPU-less, Julia-less runner; Julia/GPU-only tests are designed to skip there.CI.yml(Julia) doesn't have a GPU either, so the new CUDA-extension code path there is untested by CI (as before)🤖 Generated with Claude Code