Repository navigation
bathsearch: search every query in one pass over the target - #38
kevinmoran1988 wants to merge 11 commits into
Conversation
|
The speed/memory gains are noteworthy. This is a substantial change, so I expect it will take a while for Genevieve to work through the evaluation. I won't rush that. In the meantime, I think it's important to have a gate for very large profile sets. For example, it's not hard to imagine someone wanting to search with geNomad, which has >200,000 family profiles. |
|
Thanks, and no rush at all on the evaluation. You're right about very large profile sets. I checked how memory grows with the number of profiles, using the same 647-position profile repeated N times so that only the count changes (14 Mb target, 16 threads, peak memory):
The per-query driver stays flat and this one grows by about 0.6 MB per profile (about 1.6 MB with What I'd propose is a limit on the number of profiles searched in one pass. Profiles would be read and searched in batches of at most that many, with the target rewound between batches, and the output would stay in profile order as it is now. Two details:
If that sounds reasonable I'll add it to this PR as a separate commit. To be clear about the scope of this: the risk is only in the new driver. The other open PRs (#36, #37, TravisWheelerLab/easel#6 and TravisWheelerLab/easel#7) are speedups and memory savings in their own right and don't depend on it. The per-query column above already has them in, and it stays flat, so taking them carries no risk of running out of memory with large profile sets. Separately, I've pushed one more commit here. It fixes a slowdown I missed: my timings were all at 16 threads, and with more threads than that a search with one short profile was slower than with the per-query driver (0.59 s against 0.48 s at 24 threads, 0.70 s against 0.50 s at 32). Every change of state woke all the waiting workers, and they queued for the driver's lock while the reader waited behind them to hand over each block. The commit gives the reader its own condition variable and wakes one worker for new work. That search now takes 0.28 s at 16 threads and 0.32 s at 32. Searches with many profiles are unchanged or 2–3% faster, and the output is identical. I'm also chasing down a few more possible savings, and will post whatever holds up. 🤖 Generated with Claude Code |
|
I've added one more small commit. With 🤖 Generated with Claude Code |
|
I've added the gate as a commit here.
On the two details I raised:
So a set the size of geNomad should now need a few GB rather than well over 100. The size of the target is not a concern in the same way. The target is read in blocks, and each block's ORFs go into the same few slots, which are emptied and reused. So that part is a few MB whatever the size of the genome: it was 3 MB in a heap profile of the 352 Mb search. What does grow with the target is the hits, which are held until a batch's pass ends, and a smaller What batching doesn't lower is the floor under those numbers, which is the DP matrices of domain definition, the frameshift ones above all. A long profile against a long stretch of DNA needs a matrix of profile length × DNA length, and the largest I saw was 1.1 GB. That is the same in every driver, old or new. There are savings left there, but they would require shenanigans, so I've left it alone. The output is identical for every batch size I tried, in default, 🤖 Generated with Claude Code |
p7_pli_Frameshift() recorded which DNA window an ORF belongs to in orfsq->idx, the only thing the pipeline wrote into the ORF block it is given. Keep it in the per-call workspace instead (pli_tmp->orf_win[], starting at -1). The ORF block is then read-only to the pipeline, and an ORF needs no idx field. Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
Nothing calls it, and with the window index kept in the pipeline's workspace it was the last code that wrote an ORF's idx. Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
The pipeline reads only the residues, length and coordinates of an ORF, so it now takes easel's ESL_ORF_BLOCK in place of a block of ESL_SQ objects, about 75 bytes per ORF where an ESL_SQ took about 1 kB. bathsearch creates and empties its ORF blocks with the matching easel calls, and the seed-window scoring in p7_hmmwindow.c and the pipeline unit test pass their translated sequence as an ESL_ORF. This needs the easel commit that adds ESL_ORF_BLOCK. Output is unchanged. For one profile against a 352 Mb genome with 16 threads, peak memory falls from 133 MB to 66 MB and CPU time by 17%. Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
The clone took the source's allocation sizes, nalloc and ndomalloc, and copied that much. The domain definition's working trace keeps the size of the longest trace it has held, so every hit's copy was that large: about 2,000 states allocated for an average of 122 used, 34 kB a hit. The clone now allocates and copies what the trace holds. Hits are unchanged. This is a prerequisite for the one-pass bathsearch driver, which holds the hits of every query until the pass ends: with it, a search of 5,222 profiles against a 352 Mb genome peaks at 5.9 GB where it peaked at 28.0 GB. Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
Each step of a band folded all its vectors into one running maximum, so every max waited for the one before it. Bands of three or more vectors now keep separate maxima for even and odd vectors and join them when the band ends, as the NEON kernel already does with six. Bands of one or two vectors keep a single maximum; splitting there only adds register copies. Scores are unchanged. For models of 100 or more positions and sequences of 30 or more residues, SSV calls are 18-48% faster with AVX2 or SSE. Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
Frameshift domain definition runs om_fs5 through pli->ov5, which was only ever sized when the pipeline was created. Grow it next to ov3, so a pipeline created for a smaller model can run a larger one. Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
Two optional callbacks on P7_PIPELINE. fs_prepare() is called before om_fs3 is first read in a p7_pli_Frameshift() call, and fs5_get() supplies om_fs5 and gm_fs5 when a window reaches frameshift domain definition. A caller can then build those profiles only when a search gets that far. Both are NULL by default, which leaves the pipeline as it was; ov5 is grown and om_fs5 configured after the 5-codon pair is known. Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
bathsearch read and translated the whole target again for every query. It now reads and translates it once, and searches every query against each block. - The reader fills a small ring of block slots. A worker translates a slot once, both strands, into ORFs that are shared and read-only. - A unit of work is one query on one translated slot. A query's state (profiles, pipeline, hit list) is a lane used by one worker at a time, so nothing is copied per thread. A worker with nothing to run adds a lane to a query the others are waiting on. - DP memory belongs to the worker and is lent to the query for each unit. After a unit the worker's full matrices go back to their starting size, so no worker keeps the largest one it has met. - With --fs, a lane's frameshift profiles are built through the pipeline hooks, the first time a window needs them. - With --cpu 0, or in a build without threads, each block is read, translated and searched with every query in turn. Output is the same as the old threaded driver's at every thread count, including --cpu 0, where the old serial loop differed slightly from the threaded one. A target that can't be rewound (stdin, gzip) can now be searched with several queries, and i17-stdin.pl checks that instead of expecting a failure. 352 Mb genome, 16 threads, against the per-query driver: one profile takes the same time; 10 profiles 4.8 s to 2.9 s; 200 profiles 83 s to 47 s, and with --fs 112 s to 67 s. Peak memory for 200 profiles is 1.1 GB against 1.4 GB, and 3.3 GB against 3.9 GB with --fs. Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
The one-pass driver woke every waiting worker at each change of state. With more workers than work they all queued for the lock each time, and the reader waited behind them to hand over a block. One 115-position profile against a 352 Mb genome took 0.41 s with 16 threads, against 0.32 s with 12, and 0.70 s with 32. The reader now has its own condition variable, and new work wakes one worker; a worker that takes work wakes one more, which looks for work in its turn. The same search takes 0.28 s with 16 threads and 0.32 s with 32. Searches with many profiles are unchanged or 2-3% faster. Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
With --fs, each worker built its pair of 5-codon profiles once, for the longest query, and kept it for the whole pass. It now builds the pair for the query at hand when the query changes. Peak memory with --fs against a 352 Mb genome falls from 974 to 866 MB for 50 profiles and from 3.18 to 2.66 GB for 200. The time and the output are unchanged. Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
The one-pass driver held every query's profiles and hits for the whole pass, so its memory grew with the number of queries: about 0.5 MB each, and about 2 MB with --fs. It now reads and searches the queries --qbatch at a time (default 250), each batch in one pass over the target, and frees a batch before it reads the next. A search with no more queries than one batch is unchanged. A target that can't be read again (stdin, gzip) gets every query in its one pass, as before. For 2,000 profiles against a 352 Mb genome, peak memory falls from 3.1 to 1.7 GB and the time does not change. The output is the same for every batch size. Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
aec81e7 to
9b8e3f1
Compare
This follows up on #35. There I proposed a second program,
bathscan, and you suggested not adding a tool but lettingbathsearchpick the loop order itself. This PR goes one step further than a switch: the one-pass driver becomesbathsearch's only driver, so there is no mode to choose and no threshold to tune. That is possible because, with the changes it builds on, it takes the same time as the per-query driver for a single profile and less for two or more.What it builds on
I'm posting this now so the driver can be reviewed, but it can't be merged yet: #37 and TravisWheelerLab/easel#7 are prerequisites. The driver relies on the ORF block being read-only to the pipeline, and on hits' traces being allocated at their own length. CI here will fail until the easel change is in.
Until those are in, this branch carries their commits, and the one from #36, underneath its own. The commits that belong to this PR are the last three: growing
ov5to the model inp7_pli_Frameshift(), two optional hooks onP7_PIPELINEfor building frameshift profiles on demand, and the driver. Once the prerequisites are incorporated to your satisfaction, I'll rebase this on upstreammain, which will leave just those three.The timings below also have #36 and TravisWheelerLab/easel#6 on both sides.
What changes
bathsearchread and translated the whole target again for every query. It now reads and translates the target once, and searches every query against each block.p7_Pipeline_BATH()for that pair, and comes back for more. A slot is freed when every query has run on it.--fs, a lane's frameshift profiles are built through the pipeline hooks the first time a window needs them.--cpu 0, or in a build without threads, each block is read, translated and searched with every query in turn.After the pass, each query's hits are thresholded and written as before, including the splice step, which still runs per query.
Results
Per-query driver → one-pass driver, built from the same code apart from the driver commit. 352 Mb genome, 16 threads unless noted, AVX2. Median of 3 runs for the first six rows, of 2 for the 50-profile rows, single runs for 200 profiles.
--fs--fs523 small single-profile searches captured from a real run take 7.5–7.6 s at one thread where they took 7.7–8.1 s, and 7.2–7.3 s with
--cpu 0where they took 7.0–7.1 s.So a single profile costs the same, and two or more take 38–46% less time. Memory is lower from two profiles up, except
--fsat 50 profiles, which is about 7% higher. A single profile can use more: 30 MB more for the 647-position one, and 12 MB more at one thread.What is different for a user
--cpu 0now gives the threaded result. The old serial loop differed slightly from the threaded one. On my test set that was one E-value in the last digit in default mode, and 17 of 7,306 hits with--fs. The new driver gives the old threaded output at every thread count.i17-stdin.plexpected a failure there; it now checks that the piped target gives the same result as the file.# CPU timeandElapsedlines printed for each query now give the time since the pass began, because queries are no longer searched one after another.--fs, and a search of 5,222 profiles against this genome took 20 minutes and peaked at 5.9 GB.What I checked
--fsand--splicemodes, and--cpu 0matches 16 threads in all three.make checkpasses (easel 87/87, BATH 161/161), withi17-stdin.plupdated as above.What I haven't done
--cpu 0.This replaces a lot of
bathsearch.c, so if you'd rather keep the per-query driver available behind an option, or want this split or arranged differently, I'm happy to do that.🤖 Generated with Claude Code