Skip to content
Open
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
22 changes: 22 additions & 0 deletions ChangeLog
Original file line number Diff line number Diff line change
Expand Up @@ -14,6 +14,10 @@

## Interface and handling

- deterministic kinetic seed extension with --model=X --mode=K and configurable
local scoring; complete-energy downhill steps, atomic noLP loop/stack moves
and retained greedy traceback; equilibrium probability output is unsupported

- IntaRNAeval / --rri evaluates predefined RNA-RNA interactions (issue #184)

- compressed binary .agz accessibility caches for repeated screens (issue #245)
Expand Down Expand Up @@ -73,6 +77,24 @@ energies, restricted partition sums, and trackers.
################################################################################
################################################################################

261003 Alexander Mitrofanov
* IntaRNA/PredictorSeedExtensionKinetic, src/IntaRNA/Makefile.am :
+ greedily extend seeds on either side using complete interaction-energy
differences and deterministic thermodynamic or distance-weighted scores
+ support noLP macro-steps, GU restrictions, per-strand loop/span constraints,
explicit-seed validation and cached trajectory traceback
+ retain valid visited prefixes and select non-overlapping output from the
complete retained candidate set; reject unsupported ensemble statistics
* bin/CommandLineParsing :
+ expose --mode=K exclusively for --model=X and --kineticScore=A|B|C
* tests/PredictorSeedExtensionKinetic_test.cpp, tests/runKineticSeedExtension.sh,
tests/Makefile.am :
+ validate local move choices, macro-step energetics, structural constraints,
traceback, output filtering and CLI compatibility
* README.md, doc/kinetic-seed-extension.md, doc/Makefile.am :
+ document the council's implementation consensus, scientific scope and
replacement of unsafe loop-only pruning with exact move enumeration

261001 Alexander Mitrofanov
* IntaRNA/PredictorEvalOnly, src/IntaRNA/Makefile.am :
+ parse colon-separated hybridDB structures with sequence/index/pair validation
Expand Down
33 changes: 33 additions & 0 deletions README.md
Original file line number Diff line number Diff line change
Expand Up @@ -729,6 +729,39 @@ minimum free energy interaction.
Putative seed interactions (used by the `H` and `M` mode) can be enumerated
and studied using the `S` mode.

### Greedy kinetic seed extension

`--model=X --mode=K` grows each available seed along a deterministic greedy
path. Every step compares feasible extensions on both sides using the complete
change in interaction energy, including accessibility, terminal penalties and
dangling ends. Only strictly negative changes are accepted. With `--outNoLP`,
crossing a loop forms its closing pair and the immediately following stack as
one atomic step; a favorable stack can therefore compensate for an unfavorable
loop. Initial seeds must satisfy the selected structural constraints.

`--kineticScore` selects the local move ranking:

| Value | Score minimized for a move with gaps `s1`, `s2` |
| --- | --- |
| `A` (default) | Complete energy change |
| `B` | Complete energy change / `(1+s1+s2)` |
| `C` (C1 in the design) | Complete energy change / `(1+2*max(s1,s2))` |

Equal scores prefer the left side, then fewer unpaired bases, then smaller
`s1`. The denominators also apply to two-pair moves. All reportable visited
states, including seeds, participate in the normal energy-ranked output;
traceback preserves the actual chosen path. `--outNoGUend`, separate query and
target loop/span limits, regions, output energy/accessibility filters and
overlap settings remain applicable.

This mode is a zippering-inspired heuristic, without a calibrated time axis or
a guarantee of the global minimum. It evaluates all feasible local moves;
loop-only energetic pruning is unsafe for complete loop-plus-stack steps.
Equilibrium probability/partition-sum outputs are rejected, as are other models
and `--noSeed`. Scores B and C are optional distance preferences, not measured
kinetic rates. See the [design and implementation plan](doc/kinetic-seed-extension.md)
for the precise algorithm, scientific limitations and validation cases.


[![up](doc/figures/icon-up.28.png) back to overview](#overview)
<br /><br />
Expand Down
2 changes: 1 addition & 1 deletion doc/Makefile.am
Original file line number Diff line number Diff line change
Expand Up @@ -5,6 +5,7 @@

EXTRA_DIST = \
conda.txt \
kinetic-seed-extension.md \
doxygen.cfg \
latex-deps/adjcalc.sty \
latex-deps/adjustbox.sty \
Expand All @@ -13,4 +14,3 @@ EXTRA_DIST = \
latex-deps/tocloft.sty \
latex-deps/trimclip.sty \
latex-deps/xtab.sty

200 changes: 200 additions & 0 deletions doc/kinetic-seed-extension.md
Original file line number Diff line number Diff line change
@@ -0,0 +1,200 @@
# Deterministic kinetic seed extension

## Council decision and scientific scope

`PredictorSeedExtensionKinetic`, selected exclusively with `--model=X --mode=K`,
implements a deterministic, seed-conditioned, downhill extension heuristic.
It starts from each seed provided by the selected seed handler, compares moves
at both ends of the current duplex, and commits one move at a time. It does not
sample transition rates, simulate elapsed time, cross uphill barriers, remove
base pairs, or guarantee the global minimum free energy interaction. A
macro-step is a coarse move of this heuristic: favorable combined energy does
not establish a barrier-free physical reaction pathway.

This scope follows the distinction between gradient walks to local minima and
kinetic dynamics in the ViennaRNA ecosystem. [RNAlocmin](https://www.tbi.univie.ac.at/RNA/BHG/RNAlocmin.html)
uses gradient walks; [Kinfold and treekin](https://www.tbi.univie.ac.at/software/)
use stochastic moves or transition-rate matrices. The name "kinetic" identifies
the proposed extension strategy, not a calibrated kinetic prediction. The
choice of a single seed structure per seed start inherits the selected seed
handler's behavior; this algorithm does not enumerate every seed conformation.

The second council resolved the first review's open questions as follows.
These decisions supersede contradictory optimization and biological claims in
the supplied `PredSeedExtKinetic.md` proposal.

## State, energy and accepted moves

A state stores the complete, ordered chain of intermolecular base pairs, its
inclusive boundaries `(i1,j1,i2,j2)`, its hybridization energy `H`, and its full
interaction energy `E`. Indices use IntaRNA's internal coordinate system:
sequence 2 is reversed and prediction-range offsets are handled by the existing
wrappers. The initial `H` includes the seed's loop energies and `getE_init()`
exactly once. Every accepted step retains its actual pairs for traceback.

Evaluate every candidate with the active `InteractionEnergy` instance:

```
E_current = energy.getE(i1, j1, i2, j2, H_current)
E_next = energy.getE(i1_next, j1_next, i2_next, j2_next, H_next)
delta = E_next - E_current
```

This includes both accessibility penalties, both accessibility-weighted dangling
ends, terminal penalties, the configured temperature and parameters, and the
configured additive energy term. Changing one boundary can alter the opposite
end's dangling contribution through its accessibility weight. Consequently,
loop-plus-stack-plus-accessibility differences alone are insufficient. Infinite
states are rejected before subtraction. A constant additive term cancels in the
difference while remaining part of reported energy.

Only moves with **strictly negative full energy difference** are eligible.
Zero and uphill moves are rejected. Stop when no eligible move remains; the
finite, growing span also guarantees termination. Output energy and accessibility
thresholds filter reported states, rather than introducing additional barriers
in the trajectory. Model-inaccessible states and maximum span limits still
make a move infeasible.

For each side, enumerate all `0 <= s1 <= m1` and `0 <= s2 <= m2`, where `mk`
is the active energy model's maximum unpaired loop size for strand `k`.
`(s1,s2)=(0,0)` is a stack. Other combinations include bulges and internal loops.
The maximum total skipped length is `m1+m2`, equal to `2*m` only when the
per-strand limits are equal. The proposal's fixed 10 is not a separate limit.

A normal move adds the complementary loop-closing pair and advances each
boundary by `sk+1`. With `--outNoLP`, a move across any nonzero loop additionally
requires the immediately following outward stack. Evaluate and commit these
two pairs **atomically**, advancing each boundary by `sk+2`. Neither the
isolated closing pair nor its energy is a separately accepted or reported state.
Check both pairs, the intervening loop, and both complete strand spans before
acceptance. Each span must fit that strand's accessibility maximum length and
the requested prediction range; no single ambiguous shared `W` is introduced.

The complete initial seed must have ordered complementary pairs and finite
per-loop energies. Under `--outNoLP`, every seed pair must already have a direct
stack neighbor. Incompatible explicit seeds are skipped, including seeds with
lonely terminal pairs; this mode does not perform a preliminary seed repair.
Together with atomic macro-steps, this preserves the no-lonely-pair invariant.

Under `--outNoGUend`, both boundaries of every nonstacking extension loop must
be non-GU. Direct stacking may temporarily expose a GU outer endpoint, as in
existing IntaRNA extension recurrences. A state is reportable only if both outer
endpoints satisfy the flag. The active energy model's internal-loop GU policy
also remains authoritative. A valid earlier state remains available if descent
later stops at an unreportable GU endpoint.

## Scores and deterministic choice

`--kineticScore=A|B|C` selects the score; `C` denotes the proposal's C1:

| Option | Score minimized | Interpretation |
| --- | --- | --- |
| `A` (default) | `delta` | Steepest decrease of the actual modeled interaction energy |
| `B` | `delta/(1+s1+s2)` | Heuristic preference per total skipped length |
| `C` | `delta/(1+2*max(s1,s2))` | Heuristic preference penalizing the longer skipped strand |

A is the default because it requires no uncalibrated length-to-time assumption.
B and C are retained as explicit alternatives for exploring the supplied
proposal; their denominators are not experimentally calibrated rates or times.
The written `+1` denominator is retained even for a two-pair macro-step: it
counts one candidate move, not the number of pairs added. Scores choose the next
move only. Reported interactions remain ranked by their modeled total energy `E`.

Compare candidates from **both** sides together. Equal scores are resolved by
left before right, then smaller `s1+s2`, then smaller `s1`. Use sufficiently wide
arithmetic for score comparison, avoiding integer division truncation. Seed
iteration and output comparison follow deterministic existing IntaRNA ordering.

## Why the proposed pruning is removed

A loop-only lower bound is not a lower bound for a loop-plus-stack macro-step:
a positive loop may be rescued by the negative following stack. The initial
review reproduced a Turner2004 example at 37 degrees Celsius with a `+0.50`
kcal/mol loop and `-3.30` kcal/mol following stack, totaling `-2.80` kcal/mol.
The full delta also contains terminal and dangling changes absent from the
proposal's filters. Accessibility at the farthest candidate endpoint is not a
certified lower bound for nearer candidates; imported accessibility values need
not obey monotonicity assumptions.

Therefore the initial implementation exhaustively evaluates all feasible moves
for the current state. It uses neither the six-base-pair-type precomputed Turner
tables nor their `mdspan`/suffix-minimum representation. The suffix-minimum
operation itself is valid, but cannot repair an invalid underlying bound.
Populating an exact full-delta table first would add storage and selection work
without avoiding those energy evaluations. Omitting that cosmetic optimization
is deliberate, rather than replacing one unsafe bound with another.

There are at most `2*(m1+1)*(m2+1)` candidate shapes per state, with constant-size
local energy updates plus full boundary-energy evaluation for each. No runtime
improvement over other predictors is claimed. Any future pruning must supply a
certified lower bound for the **complete** move under the active energy model
and pass differential tests against this exhaustive implementation.

## Reporting, traceback and compatibility

The seed and every atomically committed, reportable prefix are candidates for
normal IntaRNA output. All earlier valid prefixes remain available for suboptimal
and overlap-constrained reporting. Keeping only the last state would lose valid
GU-end prefixes; keeping only the best state for each left boundary would lose
shorter candidates needed after excluding overlap with another report.

Use existing MFE ordering, output filters, seed annotations and index conversion.
Cache the actual selected pair chain for each retained candidate, including all
macro-step pairs. Traceback must recover that chain rather than run an unrelated
minimum-energy recurrence between its boundaries. Resolve duplicate boundaries
by retaining the lowest total energy, then the lexicographically smallest full
base-pair chain on an energy tie. Reduce duplicates before feeding the ordinary
optimum collector, so no stale energy can select a different cached path. Seed
annotations include only starting seeds that passed validation. Output energies
must agree with an
independent sum over the reported chain and complete boundary terms.

The cache stores full paths: its memory cost is proportional to the sum of the
retained paths' lengths, in addition to the seed handler's storage. Retaining
prefixes is therefore not a constant-memory walk. It supports reproducible
traceback and shorter alternatives for overlap-constrained output.

Only `model=X` accepts `mode=K`; seedless operation and `kineticScore` outside K
are rejected. Existing defaults and other models remain unchanged. Ordinary
energy and minimum-energy tracker output are supported. Equilibrium partition
functions and normalized equilibrium probabilities are not defined by this
selected collection of greedy trajectories. Requests needing `Zall`, including
ensemble output and probability trackers, are rejected in CLI validation, with
an API constructor guard for `needZall`. The algorithm does not manufacture an
ensemble by summing repeated prefixes from multiple seeds.

## Implementation plan and acceptance checklist

1. Add the public predictor header and implementation, derived from
`PredictorMfe`, with owned seed handler, the existing index-offset wrappers,
validated A/B/C selection and `needZall` rejection. Register both files for
building and installation.
2. Initialize and trace each handler-provided seed; validate its full chain,
seed range, energy, strand spans and active structural constraints. Record
reportable seed states.
3. Enumerate both sides and all feasible loop shapes, form complete normal or
atomic noLP moves, recompute full candidate energy, and select the strictly
downhill winner using the chosen score and the specified tie order. Repeat
until stalling or range exhaustion.
4. Preserve complete committed paths and valid prefixes; integrate ordinary MFE
output and custom candidate lookup for overlap-constrained suboptimals.
Reconstruct seed annotations without changing the retained greedy chain.
5. Wire `--model=X --mode=K` and `--kineticScore=A|B|C` into parsing, help and
factory construction. Reject incompatible model, seedless and ensemble or
probability requests with clear diagnostics. Update README and ChangeLog.
6. Add an independent tiny-sequence reference that enumerates absolute candidate
endpoints and recomputes the entire chain energy. Compare reported energy,
coordinates and traceback against it across scores, constraints and offsets.
Include targeted regressions for a positive-loop/negative-stack rescue,
strict stopping, tie order, nonmonotone accessibility, full boundary-energy
changes, invalid explicit seeds and GU-prefix retention.
7. Run focused API tests, CLI mode/flag compatibility checks, full `make tests`,
debug validation of bounds and ownership, installed standalone-header checks
where supported, and `git diff --check`. Inspect failures; never regenerate
expected outputs merely to hide a change.
8. Review the complete diff and open a pull request documenting this scientific
scope, the deliberate pruning correction, implemented behavior, validation
results and any remaining environmental limitations.

Test and build results belong in the pull request and development record; this
checklist specifies the required work and does not imply an unrun check passed.
2 changes: 2 additions & 0 deletions src/IntaRNA/Makefile.am
Original file line number Diff line number Diff line change
Expand Up @@ -72,6 +72,7 @@ libIntaRNA_a_HEADERS = \
PredictorMfe2dSeed.h \
PredictorMfe2dSeedExtension.h \
PredictorMfe2dSeedExtensionRIblast.h \
PredictorSeedExtensionKinetic.h \
PredictorMfe2dHeuristic.h \
PredictorMfe2dHeuristicSeed.h \
PredictorMfe2dHelixBlockHeuristic.h \
Expand Down Expand Up @@ -134,6 +135,7 @@ libIntaRNA_a_SOURCES = \
PredictorMfe2dSeed.cpp \
PredictorMfe2dSeedExtension.cpp \
PredictorMfe2dSeedExtensionRIblast.cpp \
PredictorSeedExtensionKinetic.cpp \
PredictorMfe2dHeuristic.cpp \
PredictorMfe2dHeuristicSeed.cpp \
PredictorMfe2dHelixBlockHeuristic.cpp \
Expand Down
Loading
Loading