From 07f931f10f21f8adfa6ed6714197eb4098b97b99 Mon Sep 17 00:00:00 2001 From: AlexanderMitrofanov Date: Wed, 21 Oct 2026 17:51:19 +0200 Subject: [PATCH 1/2] feat: select candidate sites per seed before global output filtering --- ChangeLog | 24 +++ README.md | 26 ++- doc/analysis/output-per-seed.md | 86 ++++++++ src/IntaRNA/Makefile.am | 2 + src/IntaRNA/OutputConstraint.cpp | 2 + src/IntaRNA/OutputConstraint.h | 9 +- src/IntaRNA/OutputHandlerSeedCandidates.cpp | 100 +++++++++ src/IntaRNA/OutputHandlerSeedCandidates.h | 41 ++++ src/IntaRNA/PredictorMfe.cpp | 150 ++++++++++++++ src/IntaRNA/PredictorMfe.h | 28 +++ .../PredictorMfe2dHelixBlockHeuristicSeed.cpp | 4 + src/IntaRNA/PredictorMfe2dHeuristicSeed.cpp | 4 + .../PredictorMfe2dHeuristicSeedExtension.cpp | 4 + src/IntaRNA/PredictorMfe2dSeed.cpp | 4 + src/IntaRNA/PredictorMfe2dSeedExtension.cpp | 4 + .../PredictorMfe2dSeedExtensionRIblast.cpp | 4 + ...redictorMfeEns2dHeuristicSeedExtension.cpp | 3 + .../PredictorMfeEns2dSeedExtension.cpp | 25 ++- src/IntaRNA/PredictorMfeEns2dSeedExtension.h | 8 + src/IntaRNA/PredictorMfeEnsSeedOnly.cpp | 4 + src/IntaRNA/PredictorMfeSeedOnly.cpp | 4 + src/IntaRNA/PredictorSeedExtensionKinetic.cpp | 4 + src/IntaRNA/SeededPartitionFunction.cpp | 4 + src/IntaRNA/SeededPartitionFunction.h | 5 + src/bin/CommandLineParsing.cpp | 17 +- src/bin/CommandLineParsing.h | 3 + src/bin/IntaRNA.cpp | 20 +- tests/Makefile.am | 4 +- tests/OutputHandlerSeedCandidates_test.cpp | 40 ++++ tests/PredictorSeedOracle_test.cpp | 195 ++++++++++++++++++ tests/SeededChainOracle.h | 2 + tests/SeededPartitionFunction_test.cpp | 19 +- tests/runOutputPerSeed.sh | 95 +++++++++ 33 files changed, 935 insertions(+), 9 deletions(-) create mode 100644 doc/analysis/output-per-seed.md create mode 100644 src/IntaRNA/OutputHandlerSeedCandidates.cpp create mode 100644 src/IntaRNA/OutputHandlerSeedCandidates.h create mode 100644 tests/OutputHandlerSeedCandidates_test.cpp create mode 100644 tests/runOutputPerSeed.sh diff --git a/ChangeLog b/ChangeLog index 09158edc..cb38c575 100644 --- a/ChangeLog +++ b/ChangeLog @@ -14,6 +14,10 @@ ## Interface and handling +- --outNumberPerSeed selects each seed's best candidate sites before final + deduplication, energy/overlap/count filtering; model P ranks left-most-anchor + ensembles; zero preserves ordinary prediction (#285, #288 review) + - orient probability SVGs from the upper-left corner and highlight every MFE pair, splitting seed/MFE cell borders into orange and green (PR #259 review) @@ -221,6 +225,26 @@ energies, restricted partition sums, and trackers. ################################################################################ ################################################################################ +261021 Alexander Mitrofanov + * OutputConstraint, OutputHandlerSeedCandidates, PredictorMfe, seeded predictors + + Add bounded per-seed candidate selection before global pruning, preserving + eligible energy/traceback pairs and deduplicating sites before final filters. + * Reuse each MFE mode with an admitted single-seed adapter; preserve ordinary + statistics/trackers in one separate pass. Positive quotas can cost one + prediction per seed; zero retains the existing path. + * SeededPartitionFunction, PredictorMfeEns2d[Heuristic]SeedExtension + + Rank/report left-most-anchor site ensemble energies independently of full + site partitions, Zall and actual-pair probabilities, per #288 guidance. + * CommandLineParsing, IntaRNA, README.md, doc/analysis/output-per-seed.md + + --outNumberPerSeed=[0,1000], default zero; merge seed identities across + regions/windows and reset per region combination with --outPerRegion. + * Reject positive unseeded quotas, ignore evaluation, and apply final overlap + constraints to the retained union without refilling quotas (#285). + * PredictorSeedOracle_test.cpp, OutputHandlerSeedCandidates_test.cpp, + runOutputPerSeed.sh, library/test Makefile.am + + Independent chain-enumeration, traceback, seed-only, mode-matrix, merge, + filter, CLI and unchanged-probability regressions; register/install sources. + 261020 Alexander Mitrofanov * PredictorMfe2dSeedExtension.*, PredictorMfe2dHeuristicSeedExtension.cpp - Remove extensionLoopE(), the local loop-energy table, known-value flags, diff --git a/README.md b/README.md index 5bbfd8d0..b93a8fcc 100644 --- a/README.md +++ b/README.md @@ -1766,6 +1766,30 @@ only one prediction is reported. Heuristic predictors additionally prune candidate extensions during computation and can therefore consider fewer interaction sites. +`--outNumberPerSeed=K` first retains each seed's best `K` distinct sites, then +merges their union and applies `--outNumber`, `--outOverlap`, `--outMaxE` and +`--outDeltaE`. The default `0` preserves the existing prediction path; positive +values range from 1 to 1000. With `1`, only the best candidate generated by each +seed can enter final output. Equal site boundaries are reported once, using the +lowest-energy eligible conformation. Final filtering does not refill a seed's +quota, so fewer than `outNumber` results can remain. + +MFE models rank conformations constrained to contain the generating seed. Model +P ranks and reports the ensemble energy of conformations assigned to that seed +as their left-most anchor, using the same canonical seed ownership as its +partition computation. This can differ from the full site's ensemble energy. +Ensemble totals, prediction trackers and base-pair probabilities retain their +ordinary, unlimited-seed meaning. Heuristic modes keep their own recurrence +and extension restrictions. MFE modes run a separate prediction per admitted +seed, so enabling this option can substantially increase runtime; retained +candidate storage is bounded by `K` sites per seed. + +Quotas apply across merged regions and overlapping windows of one sequence pair. +`--outPerRegion=true` resets them for each region combination, with windows of +that combination sharing a quota. Final overlap filtering of this candidate +union supports all overlap modes across regions/windows. A positive quota +requires seeds (`--noSeed` is incompatible); evaluation ignores this option. + Furthermore, it is possible to *restrict (sub)optimal enumeration* using - `--outMaxE` : maximal energy for any interaction reported @@ -1784,7 +1808,7 @@ Furthermore, it is possible to *restrict (sub)optimal enumeration* using Overlap refers to the entire interval between the outermost intermolecular base pairs, including unpaired positions inside that interval. -*Note*: non-overlapping output is heuristic, even with exact prediction +*Note*: with `--outNumberPerSeed=0`, non-overlapping output is heuristic, even with exact prediction (`--mode=M`), and increases runtime. For each left interaction boundary, only the best right extension is retained for overlap selection. If that extension overlaps an earlier result, a shorter compatible extension may already have diff --git a/doc/analysis/output-per-seed.md b/doc/analysis/output-per-seed.md new file mode 100644 index 00000000..4d0dd0b0 --- /dev/null +++ b/doc/analysis/output-per-seed.md @@ -0,0 +1,86 @@ +# Per-seed output selection (#285, review #288) + +The [maintainer's answers](https://github.com/BackofenLab/IntaRNA/pull/288#issuecomment-6319284289) +resolve the planning checkpoint in PR #288: + +- Rank each seed's own candidates before forwarding them to overall selection. +- If different seeds select different conformations at identical site boundaries, + retain the lowest-energy eligible conformation and report the site once. +- For model P, use the ensemble assigned to the current seed as its left-most + anchor. The first-seed ownership contribution is the requested objective. +- Apply quotas across merged regions/windows, or independently per region + combination with `--outPerRegion`. Deduplicate, then apply final restrictions + without refilling quotas. Reject positive quotas with `--noSeed`; ignore them + for evaluation. + +## Implementation + +`OutputConstraint::reportMaxPerSeed` is a trailing constructor argument defaulting +to zero. CLI `--outNumberPerSeed` accepts integers in [0,1000]. Zero follows the +existing ranking, traceback, ensemble and output path with no candidate storage. + +`OutputHandlerSeedCandidates` keeps at most K distinct sites for each generating +seed, keyed by its start pair in original sequence coordinates. Rediscovering a +site in an overlapping window spends no additional quota. The selected union is +deduplicated by its outer pairs; energy and conformation are retained together. +Deterministic energy/site ordering precedes final energy-window, overlap and +count filtering. No rejected candidate is replaced by the seed's next choice. +The CLI owns the collector across calls; direct predictor API users receive the +same selection within each `predict()` call. Independent regions reset quotas +only after all their windows have completed. + +MFE prediction reuses each selected mode's own implementation with a non-owning +single-seed adapter and an unrestricted intermediate output constraint. This +preserves seed-constrained conformations, traceback and the mode's heuristic +recurrences, including S/H, S/M, X/H, X/M, X/R, X/K, B/H and seed-only modes. +The adapter reuses the admitted seed table and forwards full seed annotations. +Enabling the option costs one prediction per admitted seed. When statistics or +trackers are requested, a separate ordinary prediction supplies them exactly +once; seed-specific passes never contribute. This favors explicit, verifiable +semantics over adding seed-indexed DP state to every predictor family. Retained +candidate storage is O(K times the number of seed starts), excluding the reused +prediction workspace and retained traceback lengths. Custom MFE predictor +subclasses can override `predictPerSeed` to supply a factory preserving their +recurrence/state; incompatible inherited factories raise an explicit error when +a positive quota is requested. Complete-site spans are validated before each +seed-specific rank update, including legacy DP cells reserved for extensions. +Runtime can increase +substantially on long sequences with many seeds. + +Model P consumes the existing anchor decomposition directly. The stack-seed +kernel's optional `anchorComplete` callback supplies each canonical left-most +anchor's whole hybrid partition at each complete site before full-site +aggregation. It sums noLP support states before publishing the candidate. +Legacy bulged-seed extension uses its existing anchor-specific left/right +contributions, retaining that mode's approximation and numerical limitations. +Model P reports the selected anchor ensemble's energy, including the site's +ordinary boundary/accessibility terms; it does not substitute the full site's +energy. Ensemble reports continue to represent site boundaries and annotate +the generating seed. The stack backend's canonicalization of nested explicit +seeds is unchanged: a containing seed is not an additional owner. + +`Zall`, tracker callbacks and base-pair probability matrices remain the ordinary +prediction's quantities. Neither per-seed quotas nor final overlap/count/energy +selection changes the forward/outside objective. Exact/heuristic partition flags +are forwarded with those statistics. Probability SVG structure rendering keeps +its existing site-MFE representation. + +## Validation + +Regression coverage includes independent exhaustive chain enumeration grouped +by each contained seed for exact MFE models and by the left-most seed for exact +ensembles; quotas 1, 2 and effectively unlimited; cropped/reversed coordinates; +noLP; qualifying traceback energy; seed-only prediction; duplicate sites and +windows; final overlap/energy/count selection without quota refills; and all +currently dispatched seeded CLI modes. CLI checks also cover explicit zero, +parameter files, invalid integers, unseeded/evaluation behavior, merged versus +independent regions, deterministic multithreaded windows and unchanged actual +base-pair probabilities. + +GCC 14 / ViennaRNA 2.7.2 release and debug `make tests -j2` both passed all nine +suites (123 API cases; 529,850 release / 529,851 debug assertions). The installed +pkg-config consumer and standalone compilation of the five affected public +headers passed. Out-of-tree builds used the source `src/` include path, the +generated source-side `IntaRNA.pc` expected by current build rules, and copies of +the unchanged CLI fixtures in each build's `tests/data/`. No expected prediction +outputs were regenerated. diff --git a/src/IntaRNA/Makefile.am b/src/IntaRNA/Makefile.am index 178691d6..e4100d15 100644 --- a/src/IntaRNA/Makefile.am +++ b/src/IntaRNA/Makefile.am @@ -61,6 +61,7 @@ libIntaRNA_a_HEADERS = \ OutputHandlerEnsemble.h \ OutputHandlerHub.h \ OutputHandlerInteractionList.h \ + OutputHandlerSeedCandidates.h \ OutputHandlerText.h \ PredictionTracker.h \ PredictionTrackerHub.h \ @@ -132,6 +133,7 @@ libIntaRNA_a_SOURCES = \ OutputHandlerEnsemble.cpp \ OutputHandlerCsv.cpp \ OutputHandlerInteractionList.cpp \ + OutputHandlerSeedCandidates.cpp \ OutputHandlerText.cpp \ PredictionTrackerPairMinE.cpp \ PredictionTrackerProfileMinE.cpp \ diff --git a/src/IntaRNA/OutputConstraint.cpp b/src/IntaRNA/OutputConstraint.cpp index a5babdff..9181a455 100644 --- a/src/IntaRNA/OutputConstraint.cpp +++ b/src/IntaRNA/OutputConstraint.cpp @@ -16,6 +16,7 @@ OutputConstraint::OutputConstraint( , const bool needZall , const bool needBPs , const E_type maxED + , const size_t reportMaxPerSeed ) : reportMax(reportMax) @@ -28,6 +29,7 @@ OutputConstraint::OutputConstraint( , needZall(needZall) , needBPs(needBPs) , maxED(maxED) + , reportMaxPerSeed(reportMaxPerSeed) { if(deltaE < (E_type)0.0) throw std::runtime_error("OutputConstraint(deltaE="+toString(deltaE)+") not >= 0.0"); } diff --git a/src/IntaRNA/OutputConstraint.h b/src/IntaRNA/OutputConstraint.h index 82cb8d92..c48e63d0 100644 --- a/src/IntaRNA/OutputConstraint.h +++ b/src/IntaRNA/OutputConstraint.h @@ -62,6 +62,9 @@ class OutputConstraint //! maximal ED penalty of each interacting subsequence to be considered for output const E_type maxED; + //! candidates retained per generating seed before final output selection; zero disables + const size_t reportMaxPerSeed; + public: /** @@ -78,6 +81,9 @@ class OutputConstraint * @param noGUend whether or not inter-molecular UG base pairs are allowed at interaction ends * @param needZall whether or not Zall has to be computed for output generation * @param needBPs whether or not interaction base pairs have to be traced for output generation + * @param maxED maximal accessibility penalty per interacting subsequence + * @param reportMaxPerSeed best distinct sites retained per seed before deduplication + * and final output restrictions; zero preserves ordinary prediction */ OutputConstraint( const size_t reportMax = 1 , const ReportOverlap reportOverlap = OVERLAP_BOTH @@ -88,7 +94,8 @@ class OutputConstraint , const bool noGUend = false , const bool needZall = false , const bool needBPs = true - , const E_type maxED = Accessibility::ED_UPPER_BOUND); + , const E_type maxED = Accessibility::ED_UPPER_BOUND + , const size_t reportMaxPerSeed = 0); //! destruction virtual ~OutputConstraint(); diff --git a/src/IntaRNA/OutputHandlerSeedCandidates.cpp b/src/IntaRNA/OutputHandlerSeedCandidates.cpp new file mode 100644 index 00000000..ba63a817 --- /dev/null +++ b/src/IntaRNA/OutputHandlerSeedCandidates.cpp @@ -0,0 +1,100 @@ +#include "IntaRNA/OutputHandlerSeedCandidates.h" +#include +#include +#include + +namespace IntaRNA { +namespace { +// Compare complete candidates deterministically, including degenerate traces. +bool lessCandidate(const Interaction & a, const Interaction & b) +{ + if (a < b) return true; + if (b < a) return false; + return a.basePairs < b.basePairs; +} + +bool sameSite(const Interaction & a, const Interaction & b) +{ + return a.basePairs.front() == b.basePairs.front() && a.basePairs.back() == b.basePairs.back(); +} +} // namespace + +OutputHandlerSeedCandidates::OutputHandlerSeedCandidates(const OutputConstraint & constraint) + : OutputHandler(constraint) +{ + if (!constraint.reportMaxPerSeed) + throw std::invalid_argument("per-seed selection requires a positive quota"); +} + +void OutputHandlerSeedCandidates::add(const Interaction &) +{} + +void OutputHandlerSeedCandidates::addSeed(const Interaction::BasePair & seed, const Interaction & interaction) +{ + if (interaction.isEmpty() || E_isINF(interaction.energy) || !outConstraint.reportMax) return; + // Window workers share this collector. Forward exceptions outside the block. + std::exception_ptr failure; +#if INTARNA_MULITHREADING + #pragma omp critical(intarna_omp_seedCandidates) +#endif + { + try { + auto & list = candidates[seed]; + const auto existing = std::find_if(list.begin(), list.end(), [&](const Interaction & i) { return sameSite(i,interaction); }); + const bool improves = existing != list.end() ? lessCandidate(interaction,*existing) + : list.size() < outConstraint.reportMaxPerSeed || lessCandidate(interaction,list.back()); + if (improves) { + if (existing != list.end()) list.erase(existing); + const auto position = std::lower_bound(list.begin(),list.end(),interaction,lessCandidate); + list.insert(position,interaction); + if (list.size() > outConstraint.reportMaxPerSeed) list.pop_back(); + } + } catch (...) { failure = std::current_exception(); } + } + if (failure) std::rethrow_exception(failure); +} + +void OutputHandlerSeedCandidates::report(OutputHandler & destination) const +{ + // Merge equal sites before output filtering, preserving the winning trace. + using Site = std::pair; + std::map sites; + std::vector seeds; + for (const auto & entry : candidates) seeds.push_back(entry.first); + std::sort(seeds.begin(),seeds.end()); + for (const auto & seed : seeds) for (const auto & interaction : candidates.at(seed)) { + const Site site(interaction.basePairs.front(),interaction.basePairs.back()); + auto [it,inserted] = sites.try_emplace(site,interaction); + if (!inserted && lessCandidate(interaction,it->second)) it->second = interaction; + } + std::vector ranked; + for (const auto & site : sites) ranked.push_back(site.second); + std::sort(ranked.begin(),ranked.end(),lessCandidate); + if (ranked.empty()) return; + + // Apply final restrictions once to the selected union, without quota refills. + const auto & c = outConstraint; + const auto cutoff = std::int64_t(ranked.front().energy) + c.deltaE; + std::vector reported; + for (const auto & interaction : ranked) { + if (reported.size() >= c.reportMax || interaction.energy >= c.maxE || interaction.energy > cutoff) break; + bool overlaps = false; + for (const auto * previous : reported) { + const bool target = interaction.basePairs.front().first <= previous->basePairs.back().first + && previous->basePairs.front().first <= interaction.basePairs.back().first; + const bool query = interaction.basePairs.back().second <= previous->basePairs.front().second + && previous->basePairs.back().second <= interaction.basePairs.front().second; + overlaps = overlaps || (target && (c.reportOverlap == OutputConstraint::OVERLAP_NONE || c.reportOverlap == OutputConstraint::OVERLAP_SEQ2)) + || (query && (c.reportOverlap == OutputConstraint::OVERLAP_NONE || c.reportOverlap == OutputConstraint::OVERLAP_SEQ1)); + } + if (overlaps) continue; + destination.add(interaction); + reported.push_back(&interaction); + } +} + +void OutputHandlerSeedCandidates::clear() +{ + candidates.clear(); +} +} // namespace IntaRNA diff --git a/src/IntaRNA/OutputHandlerSeedCandidates.h b/src/IntaRNA/OutputHandlerSeedCandidates.h new file mode 100644 index 00000000..cc59ae8b --- /dev/null +++ b/src/IntaRNA/OutputHandlerSeedCandidates.h @@ -0,0 +1,41 @@ +#ifndef INTARNA_OUTPUTHANDLERSEEDCANDIDATES_H_ +#define INTARNA_OUTPUTHANDLERSEEDCANDIDATES_H_ + +#include "IntaRNA/OutputHandler.h" +#include +#include + +namespace IntaRNA { + +/** Collect each seed's best distinct sites across prediction calls/windows. + * Seed keys and interactions use original sequence coordinates. Final selection + * deduplicates sites, retaining their lowest eligible energy and its structure, + * then applies energy, overlap and count restrictions without refilling quotas. + * Partition statistics are accumulated independently through OutputHandler. + */ +class OutputHandlerSeedCandidates : public OutputHandler { +public: + /** @param constraint final restrictions, including a positive per-seed quota */ + explicit OutputHandlerSeedCandidates(const OutputConstraint & constraint); + /** Ordinary ranked reports have no seed identity and are ignored. + * @param interaction ordinary report from the statistics-only prediction + */ + void add(const Interaction & interaction) override; + /** Retain a candidate before any global output pruning. + * @param seed generating seed start in original coordinates + * @param interaction candidate, including its qualifying conformation/energy + */ + void addSeed(const Interaction::BasePair & seed, const Interaction & interaction); + /** Publish the selected union, leaving partition forwarding to the caller. + * @param destination output receiving at most reportMax distinct sites + */ + void report(OutputHandler & destination) const; + /** Reset quotas for an independent region combination; retain statistics. */ + void clear(); +private: + //! Only the best K distinct sites per seed are stored, with deep-copy ownership. + std::unordered_map, Interaction::BasePair::Hash> candidates; +}; + +} // namespace IntaRNA +#endif diff --git a/src/IntaRNA/PredictorMfe.cpp b/src/IntaRNA/PredictorMfe.cpp index ac6ec7ca..b9cc1589 100644 --- a/src/IntaRNA/PredictorMfe.cpp +++ b/src/IntaRNA/PredictorMfe.cpp @@ -3,8 +3,150 @@ #include #include +#include namespace IntaRNA { +namespace { + +/** Borrow one already computed seed, while retaining all original annotations. */ +class SelectedSeed : public SeedHandlerIdxOffset { +public: + SelectedSeed(SeedHandler & original, const Interaction::BasePair & selected); + size_t fillSeed(size_t i, size_t j, size_t k, size_t l) override; + bool isSeedBound(size_t i, size_t j) const override; + E_type getSeedE(size_t i, size_t j) const override; + size_t getSeedLength1(size_t i, size_t j) const override; + size_t getSeedLength2(size_t i, size_t j) const override; + bool updateToNextSeed(size_t & i, size_t & j, size_t lo1, size_t hi1, size_t lo2, size_t hi2) const override; +private: + const Interaction::BasePair selected; +}; + +inline SelectedSeed::SelectedSeed(SeedHandler & original, const Interaction::BasePair & selected) + : SeedHandlerIdxOffset(&original,false), selected(selected) +{} + +inline bool SelectedSeed::isSeedBound(size_t i, size_t j) const +{ return selected == Interaction::BasePair(i,j) && SeedHandlerIdxOffset::isSeedBound(i,j); } +inline E_type SelectedSeed::getSeedE(size_t i, size_t j) const +{ return isSeedBound(i,j) ? SeedHandlerIdxOffset::getSeedE(i,j) : E_INF; } +inline size_t SelectedSeed::getSeedLength1(size_t i, size_t j) const +{ return isSeedBound(i,j) ? SeedHandlerIdxOffset::getSeedLength1(i,j) : 0; } +inline size_t SelectedSeed::getSeedLength2(size_t i, size_t j) const +{ return isSeedBound(i,j) ? SeedHandlerIdxOffset::getSeedLength2(i,j) : 0; } +inline size_t SelectedSeed::fillSeed(size_t i, size_t j, size_t k, size_t l) +{ return selected.first >= i && selected.first <= j && selected.second >= k && selected.second <= l && isSeedBound(selected.first,selected.second); } +inline bool SelectedSeed::updateToNextSeed(size_t & i, size_t & j, size_t lo1, size_t hi1, size_t lo2, size_t hi2) const +{ + if (Interaction::BasePair(i,j) == selected || (selected.first < lo1 || selected.first > hi1 || selected.second < lo2 || selected.second > hi2 || !isSeedBound(selected.first,selected.second))) return false; + i = selected.first; j = selected.second; + return true; +} + +/** Forward callbacks without transferring ownership of the caller's tracker. */ +class BorrowedTracker : public PredictionTracker { +public: + explicit BorrowedTracker(PredictionTracker & tracker); + void updateOptimumCalled(size_t i, size_t j, size_t k, size_t l, E_type e) override; +private: + PredictionTracker & tracker; +}; +inline BorrowedTracker::BorrowedTracker(PredictionTracker & tracker) : tracker(tracker) {} +inline void BorrowedTracker::updateOptimumCalled(size_t i,size_t j,size_t k,size_t l,E_type e) +{ tracker.updateOptimumCalled(i,j,k,l,e); } + +/** Discard ordinary ranking while collecting its unchanged partition total. */ +class StatisticsOutput : public OutputHandler { +public: + explicit StatisticsOutput(const OutputConstraint & constraint); + void add(const Interaction &) override; +}; +inline StatisticsOutput::StatisticsOutput(const OutputConstraint & constraint) : OutputHandler(constraint) {} +inline void StatisticsOutput::add(const Interaction &) {} + +/** Tag the independently ranked candidates with their generating seed. */ +class SeedOutput : public OutputHandler { +public: + SeedOutput(const OutputConstraint & constraint, OutputHandlerSeedCandidates & destination, Interaction::BasePair seed); + void add(const Interaction & interaction) override; +private: + OutputHandlerSeedCandidates & destination; + const Interaction::BasePair seed; +}; +inline SeedOutput::SeedOutput(const OutputConstraint & constraint, OutputHandlerSeedCandidates & destination, Interaction::BasePair seed) + : OutputHandler(constraint), destination(destination), seed(seed) {} +inline void SeedOutput::add(const Interaction & interaction) { destination.addSeed(seed,interaction); } + +} // namespace + +void PredictorMfe::initSeedCandidates() +{ + ownedSeedCandidates.reset(); + seedCandidates = nullptr; + if (!output.getOutputConstraint().reportMaxPerSeed) return; + seedCandidates = dynamic_cast(&output); + if (!seedCandidates) { + ownedSeedCandidates = std::make_unique(output.getOutputConstraint()); + seedCandidates = ownedSeedCandidates.get(); + } +} + +bool PredictorMfe::predictPerSeed(const IndexRange & r1, const IndexRange & r2, + SeedHandlerIdxOffset & seeds, + const std::function & factory) +{ + const auto & c = output.getOutputConstraint(); + if (!c.reportMaxPerSeed) return false; + // Validate before unsigned range arithmetic or seed-table access. + energy.setOffset1(0); energy.setOffset2(0); + const size_t n = energy.size1(), m = energy.size2(); + if (!r1.isAscending() || !r2.isAscending() || r1.from >= n || r2.from >= m) + throw std::invalid_argument("per-seed prediction: invalid sequence range"); + initSeedCandidates(); + SeedHandler & original = seeds.getOriginalSeedHandler(); + // A custom recurrence must explicitly supply its own child factory rather + // than silently reverting to a standard predictor when quotas are enabled. + auto create = [&](OutputHandler & out, PredictionTracker * tracker, SeedHandler * handler) { + std::unique_ptr child(factory(out,tracker,handler)); + if (typeid(*child) != typeid(*this)) + throw std::logic_error("per-seed prediction: custom predictor must override predictPerSeed with its own factory"); + return child; + }; + // Compute statistics once with every seed and the original prediction mode. + Zall = 0; + if (c.needZall || predTracker) { + const OutputConstraint ordinary(c.reportMax,c.reportOverlap,c.maxE,c.deltaE,c.bestSeedOnly, + c.noLP,c.noGUend,c.needZall,c.needBPs,c.maxED); + StatisticsOutput statistics(ordinary); + auto predictor = create(statistics, + predTracker ? new BorrowedTracker(*predTracker) : nullptr, + new SeedHandlerIdxOffset(&original,false)); + predictor->predict(r1,r2); + Zall = predictor->getZall(); + if (statistics.usesExactPartition()) output.setExactPartition(true); + output.incrementZ(statistics.getZ()); + } + // Reuse the admitted seed table; singleton adapters never refill or own it. + const size_t last1 = std::min(r1.to,n-1), last2 = std::min(r2.to,m-1); + if (c.reportMax && original.fillSeed(r1.from,last1,r2.from,last2)) { + const OutputConstraint perSeed(c.reportMaxPerSeed,OutputConstraint::OVERLAP_BOTH,E_INF,E_INF, + c.bestSeedOnly,c.noLP,c.noGUend,false,c.needBPs,c.maxED); + size_t i = RnaSequence::lastPos, j = RnaSequence::lastPos; + while (original.updateToNextSeed(i,j,r1.from,last1,r2.from,last2)) { + const size_t l1 = original.getSeedLength1(i,j), l2 = original.getSeedLength2(i,j); + if (!l1 || !l2 || l1-1 > last1-i || l2-1 > last2-j + || l1 > energy.getAccessibility1().getMaxLength() || l2 > energy.getAccessibility2().getMaxLength()) continue; + SeedOutput ranked(perSeed,*seedCandidates,energy.getBasePair(i,j)); + auto predictor = create(ranked,nullptr,new SelectedSeed(original,{i,j})); + // Legacy recurrences may visit cells beyond the output span. Exclude + // those before they can spend a seed's quota, retaining zero behavior. + static_cast(predictor.get())->perSeedPass = true; + predictor->predict(r1,r2); + } + } + if (ownedSeedCandidates) seedCandidates->report(output); + return true; +} //////////////////////////////////////////////////////////////////////////// @@ -34,6 +176,7 @@ void PredictorMfe:: initOptima() { + initSeedCandidates(); // temporary access const OutputConstraint & outConstraint = output.getOutputConstraint(); @@ -79,6 +222,8 @@ isValidOutputSite( const size_t i1, const size_t j1, const size_t i2, const size_t j2 ) const { const OutputConstraint & constraint = output.getOutputConstraint(); + if (perSeedPass && (j1-i1+1 > energy.getAccessibility1().getMaxLength() + || j2-i2+1 > energy.getAccessibility2().getMaxLength())) return false; if (constraint.noGUend && (energy.isGU(i1,i2) || energy.isGU(j1,j2))) { return false; } @@ -214,6 +359,11 @@ reportOptima() if (Z_isNotINF(getZall())) { output.incrementZ( getZall() ); } + // Ensemble candidates were ranked by anchor before shared site pruning. + if (seedCandidates) { + if (ownedSeedCandidates) seedCandidates->report(output); + return; + } // A zero ranked-output count still publishes ensemble statistics. if (mfeInteractions.empty()) return; diff --git a/src/IntaRNA/PredictorMfe.h b/src/IntaRNA/PredictorMfe.h index 35501d33..e272bb59 100644 --- a/src/IntaRNA/PredictorMfe.h +++ b/src/IntaRNA/PredictorMfe.h @@ -6,8 +6,12 @@ #include "IntaRNA/Predictor.h" #include "IntaRNA/IndexRangeList.h" +#include "IntaRNA/OutputHandlerSeedCandidates.h" +#include "IntaRNA/SeedHandlerIdxOffset.h" +#include #include +#include #include #include @@ -24,6 +28,30 @@ namespace IntaRNA { class PredictorMfe : public Predictor { protected: + /** Predict independently with each admitted seed, preserving this mode's + * recurrences and qualifying traceback. Statistics/trackers use one ordinary + * prediction; seed-specific passes do not contribute to them. + * @param r1 target region in original coordinates + * @param r2 reversed-query region in original coordinates + * @param seeds owned seed wrapper of this predictor + * @param factory constructs this predictor type, owning the supplied tracker + * and seed adapter; children have per-seed selection disabled + * @return true if positive-quota prediction was handled completely + * @throws std::logic_error if a custom predictor subclass does not override + * this hook to supply a factory preserving its custom recurrences + */ + virtual bool predictPerSeed(const IndexRange & r1, const IndexRange & r2, + SeedHandlerIdxOffset & seeds, + const std::function & factory); + + /** Initialize local selection, or use a caller's shared regional collector. */ + void initSeedCandidates(); + //! Active collector; null on the unchanged default path. + OutputHandlerSeedCandidates * seedCandidates = nullptr; + //! Local collector for API callers that do not provide a regional collector. + std::unique_ptr ownedSeedCandidates; + //! Enforce complete-site spans before ranking a seed-constrained child pass. + bool perSeedPass = false; /** * Describes the currently best interaction found for a left interaction diff --git a/src/IntaRNA/PredictorMfe2dHelixBlockHeuristicSeed.cpp b/src/IntaRNA/PredictorMfe2dHelixBlockHeuristicSeed.cpp index aaa6cef8..2a547a84 100644 --- a/src/IntaRNA/PredictorMfe2dHelixBlockHeuristicSeed.cpp +++ b/src/IntaRNA/PredictorMfe2dHelixBlockHeuristicSeed.cpp @@ -38,6 +38,10 @@ PredictorMfe2dHelixBlockHeuristicSeed:: predict( const IndexRange & r1 , const IndexRange & r2 ) { + if (predictPerSeed(r1,r2,seedHandler,[&](OutputHandler & out, PredictionTracker * tracker, SeedHandler * seed) { + return new PredictorMfe2dHelixBlockHeuristicSeed(energy,out,tracker,helixHandler.getConstraint(),seed); + })) return; + // temporary access const OutputConstraint & outConstraint = output.getOutputConstraint(); #if INTARNA_MULITHREADING diff --git a/src/IntaRNA/PredictorMfe2dHeuristicSeed.cpp b/src/IntaRNA/PredictorMfe2dHeuristicSeed.cpp index 2c7a29f5..f721b9bd 100644 --- a/src/IntaRNA/PredictorMfe2dHeuristicSeed.cpp +++ b/src/IntaRNA/PredictorMfe2dHeuristicSeed.cpp @@ -35,6 +35,10 @@ PredictorMfe2dHeuristicSeed:: predict( const IndexRange & r1 , const IndexRange & r2 ) { + if (predictPerSeed(r1,r2,seedHandler,[&](OutputHandler & out, PredictionTracker * tracker, SeedHandler * seed) { + return new PredictorMfe2dHeuristicSeed(energy,out,tracker,seed); + })) return; + #if INTARNA_MULITHREADING #pragma omp critical(intarna_omp_logOutput) #endif diff --git a/src/IntaRNA/PredictorMfe2dHeuristicSeedExtension.cpp b/src/IntaRNA/PredictorMfe2dHeuristicSeedExtension.cpp index df2d7bd1..3d76ea01 100644 --- a/src/IntaRNA/PredictorMfe2dHeuristicSeedExtension.cpp +++ b/src/IntaRNA/PredictorMfe2dHeuristicSeedExtension.cpp @@ -33,6 +33,10 @@ void PredictorMfe2dHeuristicSeedExtension:: predict( const IndexRange & r1, const IndexRange & r2 ) { + if (predictPerSeed(r1,r2,seedHandler,[&](OutputHandler & out, PredictionTracker * tracker, SeedHandler * seed) { + return new PredictorMfe2dHeuristicSeedExtension(energy,out,tracker,seed); + })) return; + #if INTARNA_MULITHREADING #pragma omp critical(intarna_omp_logOutput) #endif diff --git a/src/IntaRNA/PredictorMfe2dSeed.cpp b/src/IntaRNA/PredictorMfe2dSeed.cpp index acd19a75..b9d96a01 100644 --- a/src/IntaRNA/PredictorMfe2dSeed.cpp +++ b/src/IntaRNA/PredictorMfe2dSeed.cpp @@ -32,6 +32,10 @@ void PredictorMfe2dSeed:: predict( const IndexRange & r1, const IndexRange & r2 ) { + if (predictPerSeed(r1,r2,seedHandler,[&](OutputHandler & out, PredictionTracker * tracker, SeedHandler * seed) { + return new PredictorMfe2dSeed(energy,out,tracker,seed); + })) return; + #if INTARNA_MULITHREADING #pragma omp critical(intarna_omp_logOutput) #endif diff --git a/src/IntaRNA/PredictorMfe2dSeedExtension.cpp b/src/IntaRNA/PredictorMfe2dSeedExtension.cpp index c2b884c2..86b273de 100644 --- a/src/IntaRNA/PredictorMfe2dSeedExtension.cpp +++ b/src/IntaRNA/PredictorMfe2dSeedExtension.cpp @@ -33,6 +33,10 @@ void PredictorMfe2dSeedExtension:: predict( const IndexRange & r1, const IndexRange & r2 ) { + if (predictPerSeed(r1,r2,seedHandler,[&](OutputHandler & out, PredictionTracker * tracker, SeedHandler * seed) { + return new PredictorMfe2dSeedExtension(energy,out,tracker,seed); + })) return; + #if INTARNA_MULITHREADING #pragma omp critical(intarna_omp_logOutput) #endif diff --git a/src/IntaRNA/PredictorMfe2dSeedExtensionRIblast.cpp b/src/IntaRNA/PredictorMfe2dSeedExtensionRIblast.cpp index 8193906c..ae6a2741 100644 --- a/src/IntaRNA/PredictorMfe2dSeedExtensionRIblast.cpp +++ b/src/IntaRNA/PredictorMfe2dSeedExtensionRIblast.cpp @@ -32,6 +32,10 @@ void PredictorMfe2dSeedExtensionRIblast:: predict( const IndexRange & r1, const IndexRange & r2 ) { + if (predictPerSeed(r1,r2,seedHandler,[&](OutputHandler & out, PredictionTracker * tracker, SeedHandler * seed) { + return new PredictorMfe2dSeedExtensionRIblast(energy,out,tracker,seed); + })) return; + // temporary access const OutputConstraint & outConstraint = output.getOutputConstraint(); #if INTARNA_MULITHREADING diff --git a/src/IntaRNA/PredictorMfeEns2dHeuristicSeedExtension.cpp b/src/IntaRNA/PredictorMfeEns2dHeuristicSeedExtension.cpp index 84faaab0..e8fa4876 100644 --- a/src/IntaRNA/PredictorMfeEns2dHeuristicSeedExtension.cpp +++ b/src/IntaRNA/PredictorMfeEns2dHeuristicSeedExtension.cpp @@ -210,6 +210,7 @@ fillHybridZ_right( const size_t sj1, const size_t sj2 // seed only not covered due to min-val of r1,r2 if(!selectingStackSeedRight) { updateZ(si1, sj1+r1, si2, sj2+r2, seedZ * rightExtZ * initZ, true); + if (seedCandidates) addSeedPartition({si1,si2},{si1,si2},{sj1+r1,sj2+r2},seedZ*rightExtZ*initZ); addExtensionProbabilityRoot(si1,sj1+r1,si2,sj2+r2,seedZ*rightExtZ*initZ,false,true); } } @@ -246,6 +247,7 @@ fillHybridZ_left( const size_t si1, const size_t si2 ) // Z( left + seed ); covers seed only updateZ(si1-l1, sj1, si2-l2, sj2, curZ * seedZ, true); + if (seedCandidates) addSeedPartition({si1,si2},{si1-l1,si2-l2},{sj1,sj2},curZ*seedZ); addExtensionProbabilityRoot(si1-l1,sj1,si2-l2,sj2,curZ*seedZ,true,false); // Z( left + seed + rightOpt ) and rightOpt true seed extension @@ -258,6 +260,7 @@ fillHybridZ_left( const size_t si1, const size_t si2 ) ) { updateZ(si1-l1, j1opt, si2-l2, j2opt, curZ * seedZ * rightOptZ, true); + if (seedCandidates) addSeedPartition({si1,si2},{si1-l1,si2-l2},{j1opt,j2opt},curZ*seedZ*rightOptZ); addExtensionProbabilityRoot(si1-l1,j1opt,si2-l2,j2opt,curZ*seedZ*rightOptZ,true,true); } } diff --git a/src/IntaRNA/PredictorMfeEns2dSeedExtension.cpp b/src/IntaRNA/PredictorMfeEns2dSeedExtension.cpp index 33a5d99d..3f41521d 100644 --- a/src/IntaRNA/PredictorMfeEns2dSeedExtension.cpp +++ b/src/IntaRNA/PredictorMfeEns2dSeedExtension.cpp @@ -309,7 +309,9 @@ predict( const IndexRange & r1, const IndexRange & r2 ) // check complementarity of boundary if (Z_equal(hybridZ_right(r1,r2),0.0)) continue; // compute partition function given the current seed - updateZ(si1-l1, sj1+r1, si2-l2, sj2+r2, hybridZ_left(l1,l2) * seedZ * hybridZ_right(r1,r2), true); + const Z_type hybrid = hybridZ_left(l1,l2) * seedZ * hybridZ_right(r1,r2); + updateZ(si1-l1, sj1+r1, si2-l2, sj2+r2, hybrid, true); + if (seedCandidates) addSeedPartition({si1,si2},{si1-l1,si2-l2},{sj1+r1,sj2+r2},hybrid); } // r2 } // l2 @@ -376,6 +378,9 @@ PredictorMfeEns2dSeedExtension::predictStackSeeds(const IndexRange & r1,const In [&](K::Pair p,K::Pair q,Z_type h,Z_type b){updateExactCompleteZ(p.first,q.first,p.second,q.second,h,b);}, stackSeedRightSelector() }; + if (seedCandidates) weights.anchorComplete = [&](K::Pair seed,K::Pair left,K::Pair right,Z_type hybrid) { + addSeedPartition(seed,left,right,hybrid); + }; const auto result=K::compute(domain,seeds,weights,pairProbabilities!=nullptr); // Custom updateZ overrides can change the objective even in release builds. PartitionArithmetic::check(Zall); @@ -403,6 +408,24 @@ PredictorMfeEns2dSeedExtension::predictStackSeeds(const IndexRange & r1,const In } } +void PredictorMfeEns2dSeedExtension::addSeedPartition(Interaction::BasePair seed, + Interaction::BasePair left, Interaction::BasePair right, Z_type hybridZ) +{ + // Keep the mode's left-most-anchor ownership, independent of full-site Z. + if (!seedCandidates || hybridZ <= 0 || !isValidOutputSite(left.first,right.first,left.second,right.second)) return; + PartitionArithmetic::check(hybridZ); + Interaction candidate(energy.getAccessibility1().getSequence(),energy.getAccessibility2().getAccessibilityOrigin().getSequence()); + candidate.basePairs.push_back(energy.getBasePair(left.first,left.second)); + if (left != right) candidate.basePairs.push_back(energy.getBasePair(right.first,right.second)); + candidate.energy = energy.getE(left.first,right.first,left.second,right.second,energy.getE(hybridZ)); + // Ensemble output represents boundaries; annotate the anchor explicitly. + const size_t end1=seed.first+seedHandler.getSeedLength1(seed.first,seed.second)-1; + const size_t end2=seed.second+seedHandler.getSeedLength2(seed.first,seed.second)-1; + candidate.setSeedRange(energy.getBasePair(seed.first,seed.second),energy.getBasePair(end1,end2), + energy.getE(seed.first,end1,seed.second,end2,seedHandler.getSeedE(seed.first,seed.second)+energy.getE_init())); + seedCandidates->addSeed(energy.getBasePair(seed.first,seed.second),candidate); +} + E_type PredictorMfeEns2dSeedExtension:: getNonOverlappingEnergy( const size_t si1, const size_t si2, const size_t si1p, const size_t si2p ) { diff --git a/src/IntaRNA/PredictorMfeEns2dSeedExtension.h b/src/IntaRNA/PredictorMfeEns2dSeedExtension.h index 5c509295..a3e17c06 100644 --- a/src/IntaRNA/PredictorMfeEns2dSeedExtension.h +++ b/src/IntaRNA/PredictorMfeEns2dSeedExtension.h @@ -30,6 +30,14 @@ namespace IntaRNA { class PredictorMfeEns2dSeedExtension: public PredictorMfeEns { protected: + /** Retain an anchor-specific site ensemble before full-site aggregation. + * @param seed local anchor start, sequence 2 reversed + * @param left local left site boundary + * @param right local right site boundary + * @param hybridZ partition assigned to this left-most seed, including initiation + */ + void addSeedPartition(Interaction::BasePair seed, Interaction::BasePair left, + Interaction::BasePair right, Z_type hybridZ); //! matrix type to hold the partition functions for interaction site starts typedef Matrix Z2dMatrix; diff --git a/src/IntaRNA/PredictorMfeEnsSeedOnly.cpp b/src/IntaRNA/PredictorMfeEnsSeedOnly.cpp index c5e79d6d..b133bc0f 100644 --- a/src/IntaRNA/PredictorMfeEnsSeedOnly.cpp +++ b/src/IntaRNA/PredictorMfeEnsSeedOnly.cpp @@ -32,6 +32,10 @@ void PredictorMfeEnsSeedOnly:: predict( const IndexRange & r1, const IndexRange & r2 ) { + if (predictPerSeed(r1,r2,seedHandler,[&](OutputHandler & out, PredictionTracker * tracker, SeedHandler * seed) { + return new PredictorMfeEnsSeedOnly(energy,out,tracker,seed); + })) return; + #if INTARNA_MULITHREADING #pragma omp critical(intarna_omp_logOutput) #endif diff --git a/src/IntaRNA/PredictorMfeSeedOnly.cpp b/src/IntaRNA/PredictorMfeSeedOnly.cpp index afab7776..d24da3e0 100644 --- a/src/IntaRNA/PredictorMfeSeedOnly.cpp +++ b/src/IntaRNA/PredictorMfeSeedOnly.cpp @@ -32,6 +32,10 @@ void PredictorMfeSeedOnly:: predict( const IndexRange & r1, const IndexRange & r2 ) { + if (predictPerSeed(r1,r2,seedHandler,[&](OutputHandler & out, PredictionTracker * tracker, SeedHandler * seed) { + return new PredictorMfeSeedOnly(energy,out,tracker,seed); + })) return; + #if INTARNA_MULITHREADING #pragma omp critical(intarna_omp_logOutput) #endif diff --git a/src/IntaRNA/PredictorSeedExtensionKinetic.cpp b/src/IntaRNA/PredictorSeedExtensionKinetic.cpp index dba6cd3a..356f9a39 100644 --- a/src/IntaRNA/PredictorSeedExtensionKinetic.cpp +++ b/src/IntaRNA/PredictorSeedExtensionKinetic.cpp @@ -61,6 +61,10 @@ PredictorSeedExtensionKinetic::~PredictorSeedExtensionKinetic() void PredictorSeedExtensionKinetic::predict(const IndexRange & r1, const IndexRange & r2) { + if (predictPerSeed(r1,r2,seedHandler,[&](OutputHandler & out, PredictionTracker * tracker, SeedHandler * seed) { + return new PredictorSeedExtensionKinetic(energy,out,tracker,seed,score); + })) return; + const size_t size1 = energy.getAccessibility1().getSequence().size(); const size_t size2 = energy.getAccessibility2().getSequence().size(); if (!r1.isAscending() || !r2.isAscending() diff --git a/src/IntaRNA/SeededPartitionFunction.cpp b/src/IntaRNA/SeededPartitionFunction.cpp index dc4afff6..ac3c5bbb 100644 --- a/src/IntaRNA/SeededPartitionFunction.cpp +++ b/src/IntaRNA/SeededPartitionFunction.cpp @@ -237,6 +237,7 @@ computeAnchored(const Kernel::Domain & d,const StackSeedDomain & seeds, for(size_t h=supportBegin;h<2;++h) sourceLeft[h]=A::multiply(source,leftTotal(u,v)[h]); auto joinRight=[&](size_t i,size_t j) { const size_t endpoint=(u+length-1+i)*maxM+v+length-1+j; + Z_type anchorZ=0; for(size_t h=supportBegin;h<2;++h) for(size_t k=supportBegin;k<2;++k) { if(singleton && h==0 && k==0) continue; const Z_type l=leftTotal(u,v)[h],r=rightTotal(i,j)[k]; @@ -244,6 +245,7 @@ computeAnchored(const Kernel::Domain & d,const StackSeedDomain & seeds, const Z_type value=A::multiply(sourceLeft[h],r); if(!outside) { hybrid(b,endpoint)=A::add(hybrid(b,endpoint),value); + if(w.anchorComplete) anchorZ=A::add(anchorZ,value); continue; } const Z_type coefficient=boundary(b,endpoint); @@ -254,6 +256,8 @@ computeAnchored(const Kernel::Domain & d,const StackSeedDomain & seeds, const size_t frontier=singleton && i==0 && j==0?0:1; dRight(i,j)[k*2+frontier]=A::add(dRight(i,j)[k*2+frontier],A::multiply(adj,l)); } + // Publish the whole anchor-site ensemble before seeds are merged. + if(!outside && w.anchorComplete && anchorZ>0) w.anchorComplete(s,p,rightPair(i,j),anchorZ); }; if(!w.selectRight || (u==0 && v==0)) { // Exact prediction and a heuristic seed with no left extension diff --git a/src/IntaRNA/SeededPartitionFunction.h b/src/IntaRNA/SeededPartitionFunction.h index 017f1f12..42ad8439 100644 --- a/src/IntaRNA/SeededPartitionFunction.h +++ b/src/IntaRNA/SeededPartitionFunction.h @@ -44,6 +44,11 @@ class SeededPartitionFunction { * or not strictly later on both strands, contribute no two-sided extension. */ std::function selectRight; + /** Optional positive partition for one canonical left-most anchor and + * complete site (anchor, left boundary, right boundary, hybrid partition). + * Independent of complete/outside callbacks and includes initiation once. + */ + std::function anchorComplete; }; /** Raw unnormalized result; pair masses remain empty in partition-only mode. */ struct Result { Z_type z=0; Matrix mass; }; diff --git a/src/bin/CommandLineParsing.cpp b/src/bin/CommandLineParsing.cpp index 6650cde4..e06f0fe6 100644 --- a/src/bin/CommandLineParsing.cpp +++ b/src/bin/CommandLineParsing.cpp @@ -207,6 +207,7 @@ CommandLineParsing::CommandLineParsing( const Personality personality ) outPrefix2streamName(), outMode("outMode", "NDCE", 'N' ), outNumber("outNumber", 0, 1000, 1), + outNumberPerSeed("outNumberPerSeed", 0, 1000, 0), outOverlap("outOverlap", "NTQB", 'B' ), outDeltaE("outDeltaE", 0.0, 100.0, 100.0), outMaxE("outMaxE", -999.0, +999.0, 0.0), @@ -985,6 +986,11 @@ CommandLineParsing::CommandLineParsing( const Personality personality ) ->notifier(boost::bind(&CommandLineParsing::validate_numberArgument,this,outNumber,_1)) , std::string("number of (sub)optimal interactions to report" " (arg in range ["+toString(outNumber.min)+","+toString(outNumber.max)+"])").c_str()) + (outNumberPerSeed.name.c_str() + , value(&(outNumberPerSeed.val)) + ->default_value(outNumberPerSeed.def) + ->notifier(boost::bind(&CommandLineParsing::validate_numberArgument,this,outNumberPerSeed,_1)) + , "Best distinct sites retained per seed before final energy/overlap/count filtering (0 disables, range [0,1000]); model P uses left-most-anchor ensemble energies.") (outOverlap.name.c_str() , value(&(outOverlap.val)) ->default_value(outOverlap.def) @@ -1636,7 +1642,7 @@ parse(int argc, char** argv) if (windowWidth.val <= windowOverlap.val) { throw error("window-based computation: --windowWidth ("+toString(windowWidth.val)+") has to exceed --windowOverlap ("+toString(windowOverlap.val)+")"); } - if (outNumber.val > 1 && outOverlap.val != 'B') { + if (outNumber.val > 1 && outOverlap.val != 'B' && outNumberPerSeed.val == 0) { throw error("window-based computation: non-overlapping subopt output (-n > 1) only supported for --outOverlap=B"); } const bool windowNeedsZall = outMode.val == 'E' @@ -1648,6 +1654,9 @@ parse(int argc, char** argv) } + if (outNumberPerSeed.val > 0 && noSeedRequired) + throw error("--outNumberPerSeed requires seed-constrained prediction; incompatible with --noSeed"); + //////////////// MODEL //////////////////////// // final checks if parameter set compatible with model @@ -1729,7 +1738,7 @@ CommandLineParsing::prepareEvaluation( boost::program_options::variables_map & v "model", "mode", "kineticScore", "noSeed", "intLenMax", "qIntLenMax", "tIntLenMax", "intLoopMax", "qIntLoopMax", "tIntLoopMax", "qRegion", "tRegion", "qRegionLenMax", "tRegionLenMax", "windowWidth", "windowOverlap", - "outNumber", "outOverlap", "outMaxE", "outDeltaE", "outMinPu", + "outNumber", "outNumberPerSeed", "outOverlap", "outMaxE", "outDeltaE", "outMinPu", "outNoLP", "outNoGUend", "outBestSeedOnly", "outPerRegion", "outPairwise" }; bool removedSomeOptions = false; @@ -1749,6 +1758,7 @@ CommandLineParsing::prepareEvaluation( boost::program_options::variables_map & v intLenMax.val = intLenMax.def = 0; qIntLenMax.val = tIntLenMax.val = 0; windowWidth.val = 0; + outNumberPerSeed.val = 0; qRegionLenMax.val = tRegionLenMax.val = 0; outNoLP = outNoGUend = outPerRegion = outPairwise = false; // inform user about the ignored options and the retained settings @@ -1972,7 +1982,7 @@ void CommandLineParsing:: validateRegionOverlap( const IndexRangeList & ranges, const bool isQuery ) const { - if (isEvaluation() || outPerRegion || ranges.size() <= 1 || outOverlap.val == 'B') { + if (isEvaluation() || outPerRegion || outNumberPerSeed.val > 0 || ranges.size() <= 1 || outOverlap.val == 'B') { return; } // Splitting one sequence allows separate predictions to reuse the other. @@ -2261,6 +2271,7 @@ getOutputConstraint( const InteractionEnergy & energy, bool forPrediction ) con , outNeedsZall , outNeedsBPs , (outMinPu.val>0 ? std::min(Accessibility::ED_UPPER_BOUND, energy.getE( outMinPu.val )) : Accessibility::ED_UPPER_BOUND) + , outNumberPerSeed.val ); } diff --git a/src/bin/CommandLineParsing.h b/src/bin/CommandLineParsing.h index f0e796ad..a25ef17d 100644 --- a/src/bin/CommandLineParsing.h +++ b/src/bin/CommandLineParsing.h @@ -771,6 +771,9 @@ class CommandLineParsing { CharParameter outMode; //! number of (sub)optimal interactions to report NumberParameter outNumber; + + //! candidate quota per generating seed; zero disables + NumberParameter outNumberPerSeed; //! whether or not reported interactions can to be overlapping CharParameter outOverlap; //! deltaE to mfe allowed to report an interaction diff --git a/src/bin/IntaRNA.cpp b/src/bin/IntaRNA.cpp index f600484b..8ae31d62 100644 --- a/src/bin/IntaRNA.cpp +++ b/src/bin/IntaRNA.cpp @@ -23,6 +23,7 @@ INITIALIZE_EASYLOGGINGPP #include "IntaRNA/Predictor.h" #include "IntaRNA/OutputHandler.h" #include "IntaRNA/OutputHandlerInteractionList.h" +#include "IntaRNA/OutputHandlerSeedCandidates.h" using namespace IntaRNA; @@ -243,6 +244,11 @@ int main(int argc, char **argv){ (parameters.reportBestPerRegion() ? std::numeric_limits::max() : 1 ) * parameters.getOutputConstraint(*energy,true).reportMax ); + // Positive quotas retain seed identities until the agreed merge scope. + std::unique_ptr seedInteractions; + if (bestInteractions.getOutputConstraint().reportMaxPerSeed) + seedInteractions = std::make_unique(bestInteractions.getOutputConstraint()); + // One owner spans all disjoint searched regions for this sequence pair. std::unique_ptr pairProbabilities(parameters.getBasePairProbabilityResult(*energy)); @@ -292,7 +298,7 @@ int main(int argc, char **argv){ // The prediction collector merges equal-energy structures with // identical boundaries and pair counts, even if inner pairs differ. std::unique_ptr predictor(parameters.getPredictor( *energy, - parameters.isEvaluation() ? *output : bestInteractions, pairProbabilities.get() )); + parameters.isEvaluation() ? *output : (seedInteractions ? static_cast(*seedInteractions) : bestInteractions), pairProbabilities.get() )); INTARNA_CHECK_NOT_NULL(predictor.get(),"predictor initialization failed"); // run prediction for this window combination @@ -331,8 +337,20 @@ int main(int argc, char **argv){ } // if not threadAborted #endif }} // window combinations + // Windows of one region combination share quotas and final selection. + if (seedInteractions && parameters.reportBestPerRegion()) { + seedInteractions->report(bestInteractions); + seedInteractions->clear(); + } + } // target ranges } // query ranges + if (seedInteractions) { + if (!parameters.reportBestPerRegion()) seedInteractions->report(bestInteractions); + bestInteractions.setExactPartition(seedInteractions->usesExactPartition()); + bestInteractions.incrementZ(seedInteractions->getZ()); + } + bool publish=true; #if INTARNA_MULITHREADING #pragma omp flush(threadAborted) diff --git a/tests/Makefile.am b/tests/Makefile.am index 5ac04457..3f639a6b 100644 --- a/tests/Makefile.am +++ b/tests/Makefile.am @@ -22,7 +22,8 @@ dist_check_SCRIPTS = \ runIntaRNA.sh \ runIntaRNAeval.sh \ runKineticSeedExtension.sh \ - runOutputOverlap.sh + runOutputOverlap.sh \ + runOutputPerSeed.sh # the program to build check_PROGRAMS = runApiTests @@ -76,6 +77,7 @@ runApiTests_SOURCES = \ OutputStreamHandlerSortedCsv_test.cpp \ OutputHandlerHub_test.cpp \ OutputHandlerInteractionList_test.cpp \ + OutputHandlerSeedCandidates_test.cpp \ SeedHandlerExplicit_test.cpp \ SeedHandlerNoBulge_test.cpp \ SeedHandlerMfe_test.cpp \ diff --git a/tests/OutputHandlerSeedCandidates_test.cpp b/tests/OutputHandlerSeedCandidates_test.cpp new file mode 100644 index 00000000..3cd5744f --- /dev/null +++ b/tests/OutputHandlerSeedCandidates_test.cpp @@ -0,0 +1,40 @@ +#include "catch.hpp" +#include "IntaRNA/OutputHandlerSeedCandidates.h" +#include "IntaRNA/OutputHandlerInteractionList.h" +using namespace IntaRNA; + +TEST_CASE("seed quotas precede site deduplication and final output restrictions", "[OutputPerSeed]") { + #include "testEasyLoggingSetup.icc" + RnaSequence target("target","GGGGGGGGGG"),query("query","CCCCCCCCCC"); + auto candidate=[&](size_t i,size_t j,size_t q,size_t r,E_type e) { + Interaction x(target,query);x.basePairs={{i,q},{j,r}};x.energy=e;return x; + }; + const auto a=candidate(0,3,9,6,-1000),b=candidate(1,4,8,5,-900),c=candidate(6,7,3,2,-700),d=candidate(8,9,1,0,-600); + for (size_t quota : {1,2}) for (auto overlap : {OutputConstraint::OVERLAP_BOTH,OutputConstraint::OVERLAP_NONE,OutputConstraint::OVERLAP_SEQ1,OutputConstraint::OVERLAP_SEQ2}) { + CAPTURE(quota,overlap); + OutputConstraint constraint(10,overlap,E_INF,E_INF,false,false,false,false,true,Accessibility::ED_UPPER_BOUND,quota); + OutputHandlerSeedCandidates collector(constraint); + // Equal boundaries from another seed retain the eligible better trace. + auto alternative=a;alternative.energy=-950;alternative.basePairs.insert(alternative.basePairs.begin()+1,{2,7}); + collector.addSeed({0,9},b);collector.addSeed({0,9},a); + collector.addSeed({1,8},alternative);collector.addSeed({1,8},d); + // A second window rediscovers the same seed/site; it spends no new quota. + collector.addSeed({0,9},a);collector.addSeed({6,3},c); + OutputHandlerInteractionList output(constraint,10);collector.report(output); + REQUIRE(output.begin()!=output.end());REQUIRE((**output.begin()).basePairs==a.basePairs); + REQUIRE(output.reported()==(quota==2?(overlap==OutputConstraint::OVERLAP_BOTH?4:3):2)); + } + SECTION("a rejected best site does not refill its seed quota") { + OutputConstraint constraint(10,OutputConstraint::OVERLAP_NONE,E_INF,E_INF,false,false,false,false,true,Accessibility::ED_UPPER_BOUND,1); + OutputHandlerSeedCandidates collector(constraint);collector.addSeed({0,9},a); + collector.addSeed({1,8},b);collector.addSeed({1,8},d); + OutputHandlerInteractionList output(constraint,10);collector.report(output); + REQUIRE(output.reported()==1); + } + SECTION("energy window and final count are applied after the seed union") { + OutputConstraint constraint(1,OutputConstraint::OVERLAP_BOTH,-800,50,false,false,false,false,true,Accessibility::ED_UPPER_BOUND,2); + OutputHandlerSeedCandidates collector(constraint);collector.addSeed({0,9},a);collector.addSeed({6,3},c); + OutputHandlerInteractionList output(constraint,10);collector.report(output);REQUIRE(output.reported()==1); + collector.clear();OutputHandlerInteractionList empty(constraint,10);collector.report(empty);REQUIRE(empty.empty()); + } +} diff --git a/tests/PredictorSeedOracle_test.cpp b/tests/PredictorSeedOracle_test.cpp index 20b97de8..2f14acfb 100644 --- a/tests/PredictorSeedOracle_test.cpp +++ b/tests/PredictorSeedOracle_test.cpp @@ -1267,3 +1267,198 @@ TEST_CASE("native anchored heuristics preserve right-choice hooks and complete b } } } + +namespace { +using PerSeedSite = std::pair; +using PerSeedResults = std::map; + +PerSeedResults perSeedExpected(const InteractionEnergy & energy, const SeedOracleResult & oracle, + const size_t quota, const bool ensemble, const bool seedOnly = false) +{ + // Enumerate chains independently, then group MFE conformations by every + // contained seed, or ensemble weights by the first (left-most) seed only. + using Pair = Interaction::BasePair; + std::map> minima; + std::map> partitions; + for (const auto & interaction : oracle.interactions) { + if (seedOnly && interaction.chain.size() != 2) continue; + const auto & first = interaction.chain.front(), & last = interaction.chain.back(); + const PerSeedSite site(energy.getBasePair(first.i1,first.i2),energy.getBasePair(last.i1,last.i2)); + for (size_t k=0;ksecond = std::min(it->second,interaction.energy); + partitions[seed][site] += energy.getBoltzmannWeight(interaction.energy); + } + } + PerSeedResults result; + for (const auto & [seed,sites] : minima) { + std::vector ranked; + for (const auto & [site,mfe] : sites) { + Interaction candidate(energy.getAccessibility1().getSequence(),energy.getAccessibility2().getAccessibilityOrigin().getSequence()); + candidate.basePairs = {site.first,site.second}; + candidate.energy = mfe; + if (ensemble) { + const size_t i=energy.getIndex1(site.first),j=energy.getIndex1(site.second); + const size_t k=energy.getIndex2(site.first),l=energy.getIndex2(site.second); + const E_type boundary = energy.getE(i,j,k,l,E_type(0)); + candidate.energy = boundary + energy.getE(partitions.at(seed).at(site)/energy.getBoltzmannWeight(boundary)); + } + ranked.push_back(candidate); + } + std::sort(ranked.begin(),ranked.end()); + for (size_t k=0;ksecond=std::min(it->second,candidate.energy); + } + } + return result; +} + +PerSeedResults perSeedActual(const OutputHandlerInteractionList & output) +{ + PerSeedResults result; + for (const Interaction * candidate : output) { + REQUIRE(candidate->isValid()); + const auto [it,inserted]=result.emplace(PerSeedSite(candidate->basePairs.front(),candidate->basePairs.back()),candidate->energy); + REQUIRE(inserted); + } + return result; +} +} // namespace + +TEST_CASE("per-seed exact MFE and anchor ensembles match independent chain enumeration", "[PredictorSeedOracle][OutputPerSeed]") { + #include "testEasyLoggingSetup.icc" + RnaSequence target("target","AGCGGCG"),query("query","CCGCCGC"); + AccessibilityDisabled tAcc(target,6,nullptr),qAcc(query,6,nullptr); + ReverseAccessibility reverse(qAcc); + InteractionEnergyBasePair energy(tAcc,reverse,2,2,false,1,-100,3,37); + const auto sc=makeSeedConstraint(2,0,0,0,false); + for (bool noLP : {false,true}) for (bool cropped : {false,true}) for (size_t quota : {1,2,1000}) { + CAPTURE(noLP,cropped,quota); + const IndexRange r1=cropped?IndexRange(1,5):IndexRange(0,6); + const IndexRange r2=cropped?IndexRange(1,6):IndexRange(0,6); + const OutputConstraint c(1000,OutputConstraint::OVERLAP_BOTH,E_INF,E_INF,false,noLP,false,true,true,Accessibility::ED_UPPER_BOUND,quota); + const auto oracle=enumerateSeededInteractions(energy,sc,c,r1,r2); + REQUIRE_FALSE(oracle.interactions.empty()); + for (char model : {'S','X','P'}) { + CAPTURE(model); + OutputHandlerInteractionList output(c,1000); + std::unique_ptr predictor; + if (model=='S') predictor=std::make_unique(energy,output,nullptr,new SeedHandlerMfe(energy,sc)); + if (model=='X') predictor=std::make_unique(energy,output,nullptr,new SeedHandlerMfe(energy,sc)); + if (model=='P') predictor=std::make_unique(energy,output,nullptr,new SeedHandlerMfe(energy,sc)); + predictor->predict(r1,r2); + const auto actual=perSeedActual(output), expected=perSeedExpected(energy,oracle,quota,model=='P'); + std::ostringstream diagnostic; + for (const auto & [site,e] : actual) diagnostic << " actual " << site.first.first << ',' << site.first.second << '-' << site.second.first << ',' << site.second.second << ':' << e; + for (const auto & [site,e] : expected) diagnostic << " expected " << site.first.first << ',' << site.first.second << '-' << site.second.first << ',' << site.second.second << ':' << e; + CAPTURE(diagnostic.str()); + REQUIRE(actual==expected); + if (model=='P') REQUIRE(double(predictor->getZall())==Approx(double(oracle.unionPartition)).epsilon(1e-12)); + if (model!='P') for (const Interaction * interaction : output) { + // Each retained structure must be a physical chain of its stated + // energy, even when another seed has a better trace at that site. + bool found=false; + for (const auto & chain : oracle.interactions) { + Interaction::PairingVec pairs; + for (const auto & p : chain.chain) pairs.push_back(energy.getBasePair(p.i1,p.i2)); + found=found || (pairs==interaction->basePairs && chain.energy==interaction->energy); + } + REQUIRE(found); + } + } + } +} + +TEST_CASE("per-seed seed-only ranking retains every seed site", "[PredictorSeedOracle][OutputPerSeed]") { + #include "testEasyLoggingSetup.icc" + RnaSequence target("target","GGGGG"),query("query","CCCCC"); + AccessibilityDisabled tAcc(target,5,nullptr),qAcc(query,5,nullptr);ReverseAccessibility reverse(qAcc); + InteractionEnergyBasePair energy(tAcc,reverse,1,1); + const auto sc=makeSeedConstraint(2,0,0,0,false); + for (bool ensemble : {false,true}) { + const OutputConstraint c(1000,OutputConstraint::OVERLAP_BOTH,E_INF,E_INF,false,false,false,true,true,Accessibility::ED_UPPER_BOUND,1); + const auto oracle=enumerateSeededInteractions(energy,sc,c,IndexRange(0,4),IndexRange(0,4)); + OutputHandlerInteractionList output(c,1000); + std::unique_ptr predictor; + if (ensemble) predictor=std::make_unique(energy,output,nullptr,new SeedHandlerMfe(energy,sc)); + else predictor=std::make_unique(energy,output,nullptr,new SeedHandlerMfe(energy,sc)); + predictor->predict(); + REQUIRE(perSeedActual(output)==perSeedExpected(energy,oracle,1,false,true)); + } +} + +namespace { +class PerSeedHeuristicChoices : public PredictorMfeEns2dHeuristicSeedExtension { +public: + using PredictorMfeEns2dHeuristicSeedExtension::PredictorMfeEns2dHeuristicSeedExtension; + std::map choices() const; +}; +inline std::map PerSeedHeuristicChoices::choices() const +{ return {stackSeedRight.begin(),stackSeedRight.end()}; } + +PerSeedResults expectedAnchorSites(const InteractionEnergy & energy, const seeded_oracle::Result & oracle, size_t quota) +{ + PerSeedResults expected; + for (const auto & [seed,sites] : oracle.anchors) { + std::vector ranked; + for (const auto & [site,z] : sites) { + Interaction candidate(energy.getAccessibility1().getSequence(),energy.getAccessibility2().getAccessibilityOrigin().getSequence()); + candidate.basePairs={energy.getBasePair(site[0],site[1]),energy.getBasePair(site[2],site[3])}; + candidate.energy=energy.getE(Z_type(z)); // disabled accessibility and zero boundary terms + ranked.push_back(candidate); + } + std::sort(ranked.begin(),ranked.end()); + for (size_t i=0;isecond=std::min(it->second,c.energy); + } + } + return expected; +} +} + +TEST_CASE("per-seed ensemble selection uses anchor partitions in heuristic and bulged modes", "[PredictorSeedOracle][OutputPerSeed]") { + #include "testEasyLoggingSetup.icc" + using namespace seeded_oracle; + RnaSequence target("target","GGGGGG"),query("query","CCCCCC"); + AccessibilityDisabled tAcc(target,6,nullptr),qAcc(query,6,nullptr);ReverseAccessibility reverse(qAcc); + InteractionEnergyBasePair energy(tAcc,reverse,2,2); + auto weight=[&](E_type e) { return E_isINF(e)?0.L:static_cast(energy.getBoltzmannWeight(e)); }; + for (bool gapped : {false,true}) for (bool heuristic : {false,true}) for (size_t quota : {1,2,1000}) { + CAPTURE(gapped,heuristic,quota); + const SeedConstraint sc(gapped?3:2,0,0,0,E_INF,Accessibility::ED_UPPER_BOUND,E_INF, + IndexRangeList(),IndexRangeList(),gapped?"2|||&2||.|":"",false,false,false); + Model model{6,6,6,6,false,{},[](Pair){return true;}, + [&](Pair){return weight(energy.getE_init());}, + [&](Pair p,Pair q){return weight(energy.getE_interLeft(p.first,q.first,p.second,q.second));}, + [](Pair,Pair){return 1.L;}}; + if (gapped) model.seeds={{{1,1},{2,3},{3,4}}}; + else for (size_t i=0;i<5;++i) for (size_t j=0;j<5;++j) model.seeds.push_back({{i,j},{i+1,j+1}}); + OutputConstraint c(1000,OutputConstraint::OVERLAP_BOTH,E_INF,E_INF,false,false,false,true,true,Accessibility::ED_UPPER_BOUND,quota); + OutputHandlerInteractionList output(c,1000); + std::unique_ptr predictor; + SeedHandler * seeds=gapped?static_cast(new SeedHandlerExplicit(energy,sc)):new SeedHandlerMfe(energy,sc); + if (heuristic) predictor=std::make_unique(energy,output,nullptr,seeds); + else predictor=std::make_unique(energy,output,nullptr,seeds); + predictor->predict(); + if (heuristic) { + // Freeze the native right-endpoint choice, then independently enumerate + // physical chains admitted by that heuristic's exact retained domain. + if (gapped) model.selectRight=[](Pair,Pair){return Pair(4,5);}; + else { + const auto choices=static_cast(predictor.get())->choices(); + model.selectRight=[choices](Pair p,Pair){return choices.at(p);}; + } + } + const auto expected=enumerate(model); + REQUIRE(perSeedActual(output)==expectedAnchorSites(energy,expected,quota)); + REQUIRE(double(predictor->getZall())==Approx(double(expected.z)).epsilon(1e-12)); + } +} diff --git a/tests/SeededChainOracle.h b/tests/SeededChainOracle.h index d8b05523..44c0016a 100644 --- a/tests/SeededChainOracle.h +++ b/tests/SeededChainOracle.h @@ -28,6 +28,7 @@ struct Model { struct Result { long double z = 0, unseeded = 0, pairCount = 0; std::map boundaries; + std::map> anchors; std::map mass, spots; }; inline bool stack(Pair a, Pair b) { @@ -59,6 +60,7 @@ inline Result enumerate(const Model & model) { if (retained) { out.z += w; out.boundaries[{p.first, p.second, q.first, q.second}] += w; + if (w > 0) out.anchors[owner->front()][{p.first,p.second,q.first,q.second}] += w; out.pairCount += w * chain.size(); for (Pair pair : chain) out.mass[pair] += w; for (size_t i=p.first; i<=q.first; ++i) diff --git a/tests/SeededPartitionFunction_test.cpp b/tests/SeededPartitionFunction_test.cpp index 3c9da136..8d6c8c89 100644 --- a/tests/SeededPartitionFunction_test.cpp +++ b/tests/SeededPartitionFunction_test.cpp @@ -96,6 +96,13 @@ void compareKernels(const seeded_oracle::Model & m, std::map boundaries; std::map calls; std::map> coefficients; + std::map> anchors; + if(suffix) w.anchorComplete=[&](auto seed,auto p,auto q,auto h) { + REQUIRE(h>0); + const seeded_oracle::Boundary site{p.first,p.second,q.first,q.second}; + REQUIRE(anchors[seed].count(site)==0); + anchors[seed][site]=h; + }; w.complete=[&](auto p,auto q,auto h,auto b){ const seeded_oracle::Boundary boundary{p.first,p.second,q.first,q.second}; REQUIRE(h>0); @@ -108,6 +115,14 @@ void compareKernels(const seeded_oracle::Model & m, ++calls[boundary]; }; auto got=suffix?Kernel::compute(d,seeds,w,true):seededStack(d,seeds,w,true); + if(suffix) { + REQUIRE(anchors.size()==hybrids.anchors.size()); + for(const auto & [seed,sites]:hybrids.anchors) { + REQUIRE(anchors.at(seed).size()==sites.size()); + for(const auto & [site,h]:sites) + REQUIRE(double(anchors.at(seed).at(site))==Approx(double(h)).epsilon(2e-12)); + } + } REQUIRE(double(got.z)==Approx(double(oracle.z)).epsilon(2e-12)); for(const auto & [boundary,hybrid]:hybrids.boundaries) if(hybrid>0) REQUIRE(calls[boundary]==1); @@ -144,12 +159,14 @@ void compareKernels(const seeded_oracle::Model & m, // the same complete hybrid sum and boundary factor appear exactly once. const auto withPairs=boundaries; const auto withPairCoefficients=coefficients; - boundaries.clear();calls.clear();coefficients.clear(); + const auto withPairAnchors=anchors; + boundaries.clear();calls.clear();coefficients.clear();anchors.clear(); auto only=suffix?Kernel::compute(d,seeds,w,false):seededStack(d,seeds,w,false); REQUIRE(only.z==got.z); REQUIRE(only.mass.storageSize()==0); REQUIRE(boundaries==withPairs); REQUIRE(coefficients==withPairCoefficients); + REQUIRE(anchors==withPairAnchors); for(const auto & [boundary,count]:calls) REQUIRE(count==1); } } diff --git a/tests/runOutputPerSeed.sh b/tests/runOutputPerSeed.sh new file mode 100644 index 00000000..5b4df77b --- /dev/null +++ b/tests/runOutputPerSeed.sh @@ -0,0 +1,95 @@ +#!/usr/bin/env bash +# Per-seed selection across all dispatched modes and CLI merge/validation paths. +set -euo pipefail +bin="$INTARNABINPATH/src/bin/IntaRNA" +tmp=$(mktemp -d) +trap 'rm -rf "$tmp"' EXIT +python3 - "$bin" "$tmp" <<'PY' +import csv +import io +from pathlib import Path +import subprocess +import sys + +binary, tmp = sys.argv[1:] +common = [binary, '--energy=B', '--acc=N', '--seedBP=2', '--threads=1', + '--outMode=C', '--outCsvCols=start1,end1,start2,end2,E', + '--outNumber=1000', '--outMaxE=100', '--outDeltaE=100', + '--default-log-file=' + tmp + '/info.log'] + +def run(*args, ok=True): + overrides = {arg.split("=", 1)[0] for arg in args if arg.startswith("--")} + command = [arg for arg in common if arg.split("=", 1)[0] not in overrides] + list(args) + if any(arg.startswith("--outMode=") and arg != "--outMode=C" for arg in args): + command = [arg for arg in command if not arg.startswith("--outCsvCols=")] + process = subprocess.run(command, text=True, capture_output=True) + assert (process.returncode == 0) == ok, (args, process.stdout, process.stderr) + return process.stdout + +def rows(text): + return list(csv.DictReader(io.StringIO(text), delimiter=';')) + +for model, modes in [('S', 'HMS'), ('X', 'HMRKS'), ('B', 'H'), ('P', 'HMS')]: + for mode in modes: + args = ['-t', 'GGCGGG', '-q', 'CCCGCC', '--model=' + model, '--mode=' + mode, + '--helixMinBP=2', '--helixMaxBP=4'] + default = run(*args) + assert run(*args, '--outNumberPerSeed=0') == default, (model, mode, 'zero changed output') + one = rows(run(*args, '--outNumberPerSeed=1')) + two = rows(run(*args, '--outNumberPerSeed=2')) + assert one and two, (model, mode, one, two) + assert len(one) <= len(two), (model, mode) + assert len(rows(run(*args, '--outNumberPerSeed=2', '--outNumber=1'))) == 1 + assert not rows(run(*args, '--outNumberPerSeed=1', '--outNumber=0')) + # Quotas never change the ordinary partition/ensemble statistics. + if mode != 'K': + assert run(*args, '--outMode=E') == run(*args, '--outMode=E', '--outNumberPerSeed=1'), (model, mode, 'partition changed') + +args = ['-t', 'GGCGGG', '-q', 'CCCGCC', '--model=X', '--mode=M'] +for invalid in ['-1', '1001', '2147483648', 'abc']: + run(*args, '--outNumberPerSeed=' + invalid, ok=False) +run(*args, '--outNumberPerSeed=1', '--noSeed', ok=False) +assert run(*args, '--outNumberPerSeed=0', '--noSeed') +# Evaluation ignores prediction options before numeric validation. +run('-t', 'GG', '-q', 'CC', '--rri=1||&1||', '--outNumberPerSeed=-1') +parameter_file = Path(tmp) / 'parameters' +parameter_file.write_text('outNumberPerSeed=2\n') +assert run(*args, '--parameterFile=' + str(parameter_file)) == run(*args, '--outNumberPerSeed=2') +assert not rows(run('-t', 'AAAA', '-q', 'AAAA', '--outNumberPerSeed=1')) + +regions = ['-t', 'GGGGAAAAGGGG', '-q', 'CCCCAAAACCCC', '--model=X', '--mode=M', + '--tRegion=1-4,9-12', '--qRegion=1-4,9-12', '--outNumberPerSeed=1'] +for overlap in 'NTQB': + merged = rows(run(*regions, '--outOverlap=' + overlap)) + separate = rows(run(*regions, '--outOverlap=' + overlap, '--outPerRegion=true')) + assert merged and separate + if overlap == 'N': + assert len(merged) == 2 and len(separate) == 4, (merged, separate) + +# Overlapping windows share seed identities in original coordinates, including +# reversed query offsets. Stable results must not depend on worker count. +windows = ['-t', 'GGCGGGCGGGCGGGCG', '-q', 'CGCCCGCCCGCCCGCC', '--model=X', '--mode=M', + '--tIntLenMax=4', '--qIntLenMax=4', '--outNumberPerSeed=2'] +full = run(*windows) +for threads in ['1', '2']: + assert run(*windows, '--windowWidth=10', '--windowOverlap=4', '--threads=' + threads) == full + +# Prediction trackers observe the same ordinary candidates exactly once. +for model in 'SXP': + tracker = ['-t', 'GGCGGG', '-q', 'CCCGCC', '--model=' + model, '--mode=M', + '--out=tMinE:' + tmp + '/min.csv'] + run(*tracker) + expected = (Path(tmp) / 'min.csv').read_bytes() + run(*tracker, '--outNumberPerSeed=1') + assert (Path(tmp) / 'min.csv').read_bytes() == expected + +# Probability outputs retain the same raw ensemble when reporting is limited. +for mode in 'MH': + probability = ['-t', 'GGCGGG', '-q', 'CCCGCC', '--model=P', '--mode=' + mode, + '--out=bpProb:' + tmp + '/pairs.csv'] + run(*probability) + expected = (Path(tmp) / 'pairs.csv').read_bytes() + run(*probability, '--outNumberPerSeed=1') + assert (Path(tmp) / 'pairs.csv').read_bytes() == expected +print('Per-seed mode matrix, validation, region/window merges and probabilities passed') +PY From dcfe997a51dece659e72bfd933554ffcfc2204d5 Mon Sep 17 00:00:00 2001 From: AlexanderMitrofanov Date: Mon, 26 Oct 2026 18:26:08 +0100 Subject: [PATCH 2/2] feat(gui): expose per-seed output selection --- ChangeLog | 9 ++++++++- gui/docs/step3-frontend-ui.md | 5 ++++- gui/options.js | 8 +++++++- gui/protocol.js | 5 +++-- gui/tests/protocol.test.mjs | 21 +++++++++++++++++---- scripts/wasm/smoke-protocol.mjs | 2 +- scripts/wasm/test-browser.py | 13 +++++++++++++ 7 files changed, 53 insertions(+), 10 deletions(-) diff --git a/ChangeLog b/ChangeLog index 783a230e..3fca7ea7 100644 --- a/ChangeLog +++ b/ChangeLog @@ -16,7 +16,8 @@ - --outNumberPerSeed selects each seed's best candidate sites before final deduplication, energy/overlap/count filtering; model P ranks left-most-anchor - ensembles; zero preserves ordinary prediction (#285, #288 review) + ensembles; zero preserves ordinary prediction; the GUI exposes the quota next + to the global output count (#285, #288 review) - refine GUI seed and CSV control widths, filter examples by task, cap result table height and use dark active tabs for additional output views (PR #303 review) @@ -314,6 +315,12 @@ energies, restricted partition sums, and trackers. regions/windows and reset per region combination with --outPerRegion. * Reject positive unseeded quotas, ignore evaluation, and apply final overlap constraints to the retained union without refilling quotas (#285). + * gui/{options,protocol}.js, gui/tests/protocol.test.mjs, + scripts/wasm/{smoke-protocol.mjs,test-browser.py} + + Expose outNumberPerSeed in the first GUI output row between outNumber and + outOverlap; move energy/accessibility filters to a dedicated second row. + + Match native quota validation, window/region compatibility and JSON replay; + cover the default, range, rendered order and real WASM execution. * PredictorSeedOracle_test.cpp, OutputHandlerSeedCandidates_test.cpp, runOutputPerSeed.sh, library/test Makefile.am + Independent chain-enumeration, traceback, seed-only, mode-matrix, merge, diff --git a/gui/docs/step3-frontend-ui.md b/gui/docs/step3-frontend-ui.md index 3e41faf5..4274ccb2 100644 --- a/gui/docs/step3-frontend-ui.md +++ b/gui/docs/step3-frontend-ui.md @@ -20,7 +20,10 @@ keep their scientific grouping without stretching across the complete panel, use 16 px horizontal gaps and bottom-align labels above their controls. The Seed tab places the constraint selector (12 rem) beside Explicit seeds on desktop; the CSV sort selector uses a 19 rem track. Both controls stack at mobile -widths. +widths. Output places the global count, per-seed candidate count and overlap mode +in its first row; energy/accessibility thresholds form the second row and the +existing boolean filters follow in the third. Per-seed selection defaults to zero +and accepts the native 0–1,000 range. The institutional link reads “Bioinformatics - University of Freiburg” and opens the group website. Footer and table notes use the available width, wrapping only when the viewport requires it. diff --git a/gui/options.js b/gui/options.js index 8bc600ff..23010c7b 100644 --- a/gui/options.js +++ b/gui/options.js @@ -51,6 +51,7 @@ const fields = new Map([ number('windowWidth', 'Prediction window width', 0, 0, 99999, { ...prediction, hint: '0 disables prediction windows; otherwise at least 10 nt.' }), number('windowOverlap', 'Prediction window overlap', 150, 0, 99999, prediction), number('outNumber', 'Maximum interactions per sequence pair', 1, 0, 100, { ...prediction, hint: 'Browser limit: 100; 0 computes accessibility only.' }), + number('outNumberPerSeed', 'Best sites retained per seed', 0, 0, 1000, { ...prediction, hint: '0 disables per-seed candidate selection; positive values can increase runtime.' }), select('outOverlap', 'Allow overlapping interactions', 'B', { B: 'Both strands', N: 'Neither strand', T: 'Target only', Q: 'Query only' }, prediction), number('outMaxE', 'Maximum interaction energy (exclusive, kcal/mol)', 0, -999, 999, { ...prediction, step: 'any' }), number('outDeltaE', 'Maximum difference from MFE (kcal/mol)', 100, 0, 100, { ...prediction, step: 'any' }), @@ -86,7 +87,12 @@ export const GROUPS = [ ['noSeed'], ['seedTQ'], ['seedBP', 'seedMaxUP', 'seedQMaxUP', 'seedTMaxUP'], ['seedMaxE', 'seedMinPu', 'seedMaxEhybrid'], ['seedQRange', 'seedTRange'], ['seedNoGU', 'seedNoGUend'], ]) }, - { key: 'output-options', label: 'Output', fields: [...fields.values()].filter((field) => field.key.startsWith('out') && field.key !== 'outCsvSort') }, + { key: 'output-options', label: 'Output', rows: rows([ + ['outNumber', 'outNumberPerSeed', 'outOverlap'], + ['outMaxE', 'outDeltaE', 'outMinPu'], + ['outBestSeedOnly', 'outNoLP', 'outNoGUend', 'outPerRegion', 'outPairwise'], + ['outSpotProbSpots'], + ]) }, { key: 'csv-options', label: 'CSV', fields: [fields.get('outCsvSort')] }, { key: 'interaction-options', label: 'Model & Energy', rows: rows([ ['model', 'mode'], ['energyVRNA', 'temperature', 'energyAdd'], ['energyNoDangles'], diff --git a/gui/protocol.js b/gui/protocol.js index d452de28..d4c71c2d 100644 --- a/gui/protocol.js +++ b/gui/protocol.js @@ -187,6 +187,7 @@ export function prepareRequest(request) { } else { if (parameters.model === 'B' && parameters.mode !== 'H') throw new Error('Helix-block model B supports heuristic mode only.'); if (parameters.noSeed && ['S', 'K'].includes(parameters.mode)) throw new Error('Seed-only and greedy prediction require a seed.'); + if (parameters.noSeed && parameters.outNumberPerSeed > 0) throw new Error('Per-seed candidate selection requires a seed.'); if (parameters.mode === 'K' && (parameters.model !== 'X' || needsEnsemble)) throw new Error('Greedy prediction requires model X and cannot produce ensemble/probability output.'); if (parameters.seedTQ) validateInteractions(parameters.seedTQ, oneRecord(qRecords, 'Explicit seeds'), oneRecord(tRecords, 'Explicit seeds'), parameters.qIdxPos0, parameters.tIdxPos0, ','); const minimumSeedLength = parameters.seedTQ ? Math.min(...parameters.seedTQ.split(',').map((seed) => seed.split('&')[0].split('|').length - 1)) : parameters.seedBP; @@ -196,11 +197,11 @@ export function prepareRequest(request) { if (outputs.some((key) => ['bpProb', 'bpsvg'].includes(key)) && (parameters.model !== 'P' || !['H', 'M'].includes(parameters.mode) || parameters.windowWidth || parameters.mode === 'M' && !parameters.noSeed && (parameters.seedTQ || parameters.seedMaxUP || parameters.seedQMaxUP > 0 || parameters.seedTMaxUP > 0))) throw new Error('Base-pair probabilities require model P, heuristic/exact mode and no windows; exact mode supports no seed or stack-only seeds.'); if (parameters.windowWidth) { const maxLength = Math.max(...['q', 't'].map((p) => parameters[`${p}Acc`] === 'N' ? parameters[`${p}IntLenMax`] : Math.max(parameters[`${p}AccW`], parameters[`${p}IntLenMax`]))); - if (parameters.windowWidth < 10 || parameters.windowOverlap >= parameters.windowWidth || parameters.windowOverlap < maxLength || ['q', 't'].some((p) => !parameters[`${p}IntLenMax`] && (parameters[`${p}Acc`] === 'N' || !parameters[`${p}AccW`])) || parameters.outNumber > 1 && parameters.outOverlap !== 'B' || needsEnsemble || parameters.model === 'P') throw new Error('Prediction windows require width ≥ 10, overlap below width and at least the bounded interaction lengths, overlapping output, and no ensemble statistics/model P.'); + if (parameters.windowWidth < 10 || parameters.windowOverlap >= parameters.windowWidth || parameters.windowOverlap < maxLength || ['q', 't'].some((p) => !parameters[`${p}IntLenMax`] && (parameters[`${p}Acc`] === 'N' || !parameters[`${p}AccW`])) || parameters.outNumberPerSeed === 0 && parameters.outNumber > 1 && parameters.outOverlap !== 'B' || needsEnsemble || parameters.model === 'P') throw new Error('Prediction windows require width ≥ 10, overlap below width and at least the bounded interaction lengths, overlapping output unless per-seed selection is enabled, and no ensemble statistics/model P.'); } // Non-overlap constraints in native merged regions require a single region on // the relevant strand; automatic region detection is only safe per region. - if (!parameters.outPerRegion && parameters.outOverlap !== 'B') { + if (!parameters.outPerRegion && parameters.outNumberPerSeed === 0 && parameters.outOverlap !== 'B') { for (const p of ['q', 't']) { if ((parameters.outOverlap === 'N' || parameters.outOverlap === (p === 'q' ? 'T' : 'Q')) && (parameters[`${p}Region`].includes(',') || parameters[`${p}RegionLenMax`])) throw new Error('Non-overlapping output with multiple/automatic regions requires output per region.'); } diff --git a/gui/tests/protocol.test.mjs b/gui/tests/protocol.test.mjs index 07a60071..4fa301e6 100644 --- a/gui/tests/protocol.test.mjs +++ b/gui/tests/protocol.test.mjs @@ -2,7 +2,7 @@ import { test } from 'node:test'; import assert from 'node:assert/strict'; import { readFileSync } from 'node:fs'; import { DEFAULTS, EVALUATION_EXAMPLE, parseSequence, parseFasta, prepareRequest, makeArguments, CSV_COLUMNS, sequenceOffset } from '../protocol.js'; -import { FIELDS } from '../options.js'; +import { FIELDS, GROUPS } from '../options.js'; import { CSV_CATALOG, CSV_GROUPS, defaultColumns } from '../columns.js'; const request = (parameters = {}, extra = {}) => ({ ...DEFAULTS, parameters, ...extra }); @@ -22,12 +22,13 @@ test('FASTA normalization, identifiers, bounded multi-record work and subsets', test('CLI defaults, ranges and strict typed controls reject invalid values', () => { const prepared = prepareRequest(DEFAULTS); assert.equal(prepared.parameters.outNumber, 1); + assert.equal(prepared.parameters.outNumberPerSeed, 0); assert.equal(prepared.parameters.qIntLoopMax, 16); assert.equal(prepared.parameters.tIntLoopMax, 10); assert.equal(prepared.csvColumns.join(','), CSV_COLUMNS); assert.ok(!makeArguments(prepared).some((arg) => arg.startsWith('--personality='))); assert.ok(!makeArguments(prepared).some((arg) => arg.startsWith('--kineticScore='))); - for (const parameters of [{ outNumber: 101 }, { outNumber: 1.5 }, { temperature: '' }, { temperature: Infinity }, { qAccW: 2 }, { qAccW: 50 }, { model: '--help' }, { noSeed: 'false' }, { seedBP: 11 }]) assert.throws(() => prepareRequest(request(parameters))); + for (const parameters of [{ outNumber: 101 }, { outNumber: 1.5 }, { outNumberPerSeed: -1 }, { outNumberPerSeed: 1001 }, { outNumberPerSeed: 1.5 }, { temperature: '' }, { temperature: Infinity }, { qAccW: 2 }, { qAccW: 50 }, { model: '--help' }, { noSeed: 'false' }, { seedBP: 11 }]) assert.throws(() => prepareRequest(request(parameters))); for (const qIdxPos0 of ['--help', NaN, 2e9+1, '', false]) assert.throws(() => prepareRequest(request({}, { qIdxPos0 }))); assert.doesNotThrow(() => prepareRequest(request({}, { qIdxPos0: -2e9 }))); assert.doesNotThrow(() => prepareRequest(request({}, { qIdxPos0: 0 }))); @@ -49,8 +50,10 @@ test('evaluation validates antiparallel supplied structures in signed coordinate }); test('cross-option constraints, regions, probabilities, SHAPE and local file arguments', () => { - for (const parameters of [{ noSeed: true, mode: 'S' }, { mode: 'K', model: 'P' }, { model: 'B', mode: 'M' }, { model: 'B', helixMaxBP: 4 }, { qRegion: '8-2' }, { qAccConstr: 'b:1-20' }, { qAccConstr: 'bbb' }, { qShapeMethod: 'X' }, { tShapeConversion: 'bad' }, { qShape: '1 nope' }, { qAcc: 'P' }, { qAcc: 'P', qAccFile: 'bad' }, { outSpotProbSpots: '20&1' }, { windowWidth: 100 }, { qIntLenMax: 3 }, { qAcc: 'N', qRegionLenMax: 20 }, { qRegion: '1-10', qRegionLenMax: 20 }, { seedTQ: '1.|||||||||&1|||||||||.' }]) assert.throws(() => prepareRequest(request(parameters))); + for (const parameters of [{ noSeed: true, mode: 'S' }, { noSeed: true, outNumberPerSeed: 1 }, { mode: 'K', model: 'P' }, { model: 'B', mode: 'M' }, { model: 'B', helixMaxBP: 4 }, { qRegion: '8-2' }, { qAccConstr: 'b:1-20' }, { qAccConstr: 'bbb' }, { qShapeMethod: 'X' }, { tShapeConversion: 'bad' }, { qShape: '1 nope' }, { qAcc: 'P' }, { qAcc: 'P', qAccFile: 'bad' }, { outSpotProbSpots: '20&1' }, { windowWidth: 100 }, { qIntLenMax: 3 }, { qAcc: 'N', qRegionLenMax: 20 }, { qRegion: '1-10', qRegionLenMax: 20 }, { seedTQ: '1.|||||||||&1|||||||||.' }]) assert.throws(() => prepareRequest(request(parameters))); assert.doesNotThrow(() => prepareRequest(request({ noSeed: true, qAccConstr: 'b:1-2,p:8-9', qRegion: '1-10' }))); + assert.doesNotThrow(() => prepareRequest(request({ outNumberPerSeed: 1, outNumber: 3, outOverlap: 'N', windowWidth: 200 }))); + assert.doesNotThrow(() => prepareRequest(request({ outNumberPerSeed: 1, outOverlap: 'N', qRegion: '1-4,6-10' }))); assert.throws(() => prepareRequest(request({ mode: 'K' }, { csvColumns: ['Eall'] })), /ensemble/); assert.throws(() => prepareRequest(request({}, { outputs: ['bpsvg'] })), /model P/); assert.throws(() => prepareRequest(request({ outCsvSort: 'ED1' })), /sort column/); @@ -62,6 +65,15 @@ test('cross-option constraints, regions, probabilities, SHAPE and local file arg assert.ok(!makeArguments(prepared, 'D').some((arg) => arg.startsWith('--outCsvCols'))); }); +test('output controls preserve the reviewed three-row filter order', () => { + const rows = GROUPS.find((group) => group.key === 'output-options').rows.map((row) => row.map((field) => field.key)); + assert.deepEqual(rows.slice(0, 3), [ + ['outNumber', 'outNumberPerSeed', 'outOverlap'], + ['outMaxE', 'outDeltaE', 'outMinPu'], + ['outBestSeedOnly', 'outNoLP', 'outNoGUend', 'outPerRegion', 'outPairwise'], + ]); +}); + test('browser catalog covers every native CSV column and every control names a CLI option', () => { const header = readFileSync(new URL('../../src/IntaRNA/OutputHandlerCsv.h', import.meta.url), 'utf8'); const ids = [...header.matchAll(/colType2string\[\w+\] = "(\w+)"/g)].map((match) => match[1]); @@ -111,9 +123,10 @@ test('CSV groups cover native columns once and task defaults serialize explicitl for (const key of ['seedStart1', 'seedStart2', 'seedEnd1', 'seedEnd2']) assert.ok(!defaultColumns('evaluation').includes(key), key); for (const key of ['ED1', 'ED2', 'seedE', 'seedStart1', 'seedStart2']) assert.ok(defaultColumns('prediction').includes(key), key); const { exportParameters, importParameters } = await import('../parameter-file.js'); - const input = request({ outNumber: 10 }, { csvColumns: defaultColumns('prediction') }); + const input = request({ outNumber: 10, outNumberPerSeed: 2 }, { csvColumns: defaultColumns('prediction') }); const saved = exportParameters(input, 'test'); assert.equal(saved.parameters.outNumber, 10); + assert.equal(saved.parameters.outNumberPerSeed, 2); assert.equal(saved.parameters.outCsvCols, input.csvColumns.join(',')); assert.deepEqual(prepareRequest(importParameters(saved)), prepareRequest(input)); }); diff --git a/scripts/wasm/smoke-protocol.mjs b/scripts/wasm/smoke-protocol.mjs index 9dd90369..09c8b400 100644 --- a/scripts/wasm/smoke-protocol.mjs +++ b/scripts/wasm/smoke-protocol.mjs @@ -44,7 +44,7 @@ assert.ok(probabilities.fs.readFile('/output/bpsvg.svg', { encoding: 'utf8' }).i const accessibility = await run({ parameters: { outNumber: 0 }, outputs: ['qPu'] }); const imported = await run({ parameters: { qAcc: 'P', qAccFile: accessibility.fs.readFile('/output/qPu.csv', { encoding: 'utf8' }) } }); assert.equal(imported.rows[0].E, (await run({ parameters: { qIntLenMax: 9 } })).rows[0].E); -for (const parameters of [{ energyVRNA: 'Turner99', temperature: 20 }, { energyVRNA: 'Andronescu07' }, { model: 'B' }, { mode: 'K' }, { noSeed: true, mode: 'M' }, { seedMaxUP: 1 }, { qAccConstr: '..........', qShape: '1 0.2\n2 0.3' }, { outOverlap: 'N', outNumber: 3 }, { outCsvSort: 'E', outNumber: 3 }]) assert.ok((await run({ parameters })).rows.length >= 1); +for (const parameters of [{ energyVRNA: 'Turner99', temperature: 20 }, { energyVRNA: 'Andronescu07' }, { model: 'B' }, { mode: 'K' }, { noSeed: true, mode: 'M' }, { seedMaxUP: 1 }, { qAccConstr: '..........', qShape: '1 0.2\n2 0.3' }, { outNumberPerSeed: 1 }, { outOverlap: 'N', outNumber: 3 }, { outCsvSort: 'E', outNumber: 3 }]) assert.ok((await run({ parameters })).rows.length >= 1); console.log('Browser protocol: native-reference prediction/evaluation, signed coordinates, multi-FASTA, all columns, imported accessibility, energy/seed/helix/SHAPE options and output files passed'); // Execute exported parameter files in fresh engines with only archived inputs. diff --git a/scripts/wasm/test-browser.py b/scripts/wasm/test-browser.py index 5152f68f..cc62414a 100644 --- a/scripts/wasm/test-browser.py +++ b/scripts/wasm/test-browser.py @@ -113,6 +113,7 @@ def predict(): expect(page.locator('.brand')).to_contain_text('IntaRNA package') expect(page.locator('#results-title')).to_have_text('Results') expect(page.locator('#outNumber')).to_have_value('10') + expect(page.locator('#outNumberPerSeed')).to_have_value('0') assert page.locator('.compute-bar').evaluate('node => node.firstElementChild.id') == 'compute' assert page.locator('.sequence-meta').first.evaluate('node => node.firstElementChild.tagName') == 'LABEL' fill('query', '>named\tignored description\n' + 'G' * 1000 + '\n' * 20 + 'GG') @@ -174,6 +175,18 @@ def predict(): expect(page.locator('#seedBP')).to_be_enabled() expect(page.locator('#seedTQ')).to_be_disabled() open_group('output-options') + output_rows = [ + [page.locator('#' + key).locator('xpath=..').bounding_box() for key in keys] + for keys in [ + ['outNumber', 'outNumberPerSeed', 'outOverlap'], + ['outMaxE', 'outDeltaE', 'outMinPu'], + ['outBestSeedOnly', 'outNoLP', 'outNoGUend', 'outPerRegion', 'outPairwise'], + ] + ] + for boxes in output_rows: + assert max(box['y'] for box in boxes) - min(box['y'] for box in boxes) < 1, boxes + assert output_rows[0][0]['x'] < output_rows[0][1]['x'] < output_rows[0][2]['x'] + assert output_rows[0][0]['y'] < output_rows[1][0]['y'] < output_rows[2][0]['y'] output_labels = page.locator('#output-options .checkbox-grid label') assert output_labels.count() == 15 assert all(label.get_attribute('title') for label in output_labels.all())