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
74 changes: 67 additions & 7 deletions src/IntaRNA/PredictorMfeEns.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -14,6 +14,7 @@ PredictorMfeEns::PredictorMfeEns(
, PredictionTracker * predTracker
)
: PredictorMfe(energy,output,predTracker)
, updateZisComplete(false)
{
}

Expand All @@ -36,34 +37,35 @@ initZ()

////////////////////////////////////////////////////////////////////////////

void
bool
PredictorMfeEns::
updateZ( const size_t i1, const size_t j1
addPartitionContribution( const size_t i1, const size_t j1
, const size_t i2, const size_t j2
, const Z_type partZ
, const bool isHybridZ )
, const bool isHybridZ
, Z_type & partZ_noED )
{
// check if something to be done
if (Z_equal(partZ,0) || Z_isINF(Zall))
return;
return false;

// Apply the same site filters used for MFE candidates before changing
// either the global or boundary-specific partition.
const OutputConstraint & outConstraint = output.getOutputConstraint();
if (outConstraint.noGUend
&& (energy.isGU(i1,i2) || energy.isGU(j1,j2)))
{
return;
return false;
}
if (outConstraint.maxED < Accessibility::ED_UPPER_BOUND
&& (energy.getED1(i1,j1) > outConstraint.maxED
|| energy.getED2(i2,j2) > outConstraint.maxED))
{
return;
return false;
}

// handle whether or not partZ includes ED values or not
Z_type partZ_withED = 0, partZ_noED = 0;
Z_type partZ_withED = 0;
if (isHybridZ) {
#if INTARNA_IN_DEBUG_MODE
if ( (std::numeric_limits<Z_type>::max() - (partZ*energy.getBoltzmannWeight(energy.getE(i1,j1,i2,j2, E_type(0))))) <= Zall) {
Expand All @@ -86,6 +88,33 @@ updateZ( const size_t i1, const size_t j1

// increase overall partition function
Zall += partZ_withED;
return true;
}

////////////////////////////////////////////////////////////////////////////

void
PredictorMfeEns::
updateZ( const size_t i1, const size_t j1
, const size_t i2, const size_t j2
, const Z_type partZ
, const bool isHybridZ )
{
Z_type partZ_noED = 0;
if (!addPartitionContribution(i1, j1, i2, j2, partZ, isHybridZ, partZ_noED)) {
return;
}

// Exact 2D prediction finalizes each boundary once. In this scoped mode the
// partition can update the optima immediately instead of entering the map.
// Trackers retain the map path to preserve deferred callback ordering.
if (updateZisComplete && predTracker == NULL) {
if (Z_isNotINF(partZ_noED) && partZ_noED > 0) {
PredictorMfe::updateOptima(i1, j1, i2, j2,
energy.getE(partZ_noED), true, false);
}
return;
}

// store partial Z (without ED)
Interaction::Boundary key(i1,j1,i2,j2);
Expand All @@ -101,6 +130,37 @@ updateZ( const size_t i1, const size_t j1

////////////////////////////////////////////////////////////////////////////

void
PredictorMfeEns::
updateCompleteZ( const size_t i1, const size_t j1
, const size_t i2, const size_t j2
, const Z_type partZ
, const bool isHybridZ )
{
class CompleteUpdateScope {
public:
explicit CompleteUpdateScope(bool & mode)
: mode(mode), previous(mode)
{
mode = true;
}

~CompleteUpdateScope()
{
mode = previous;
}

private:
bool & mode;
const bool previous;
};

CompleteUpdateScope completeUpdate(updateZisComplete);
updateZ(i1, j1, i2, j2, partZ, isHybridZ);
}

////////////////////////////////////////////////////////////////////////////

void
PredictorMfeEns::
updateOptimaUsingZ()
Expand Down
32 changes: 31 additions & 1 deletion src/IntaRNA/PredictorMfeEns.h
Original file line number Diff line number Diff line change
Expand Up @@ -61,6 +61,9 @@ class PredictorMfeEns : public PredictorMfe {
//! map storing the partition of Zall for all considered interaction sites
std::unordered_map<Interaction::Boundary, Z_type, Interaction::Boundary::Hash> Z_partition;

//! whether updateZ() currently receives one complete boundary partition
bool updateZisComplete;


/**
* Initializes the hybridization partition functions.
Expand Down Expand Up @@ -94,6 +97,20 @@ class PredictorMfeEns : public PredictorMfe {
, const Z_type partFunct
, const bool isHybridZ );

/**
* Consumes a complete partition for one interaction boundary. Exact 2D
* prediction finalizes each boundary once and therefore does not need to
* retain it until reportOptima(). The virtual updateZ() hook is called so
* subclasses can observe or alter the update; overrides that delegate to the
* base implementation inherit streaming. Tracker-enabled predictions retain
* the partition so their established callback timing and ordering stay intact.
*/
void
updateCompleteZ( const size_t i1, const size_t j1
, const size_t i2, const size_t j2
, const Z_type partFunct
, const bool isHybridZ );

/**
* Calls for the stored Z_partition information updateOptima() before
* calling reportOptima() from its super class.
Expand All @@ -113,7 +130,20 @@ class PredictorMfeEns : public PredictorMfe {
private:

/**
* Calls updateOptima() for each entry of Z_partition.
* Applies output filters, converts the partition representation and adds
* the accepted contribution to Zall. Returns the partition without ED.
*/
bool
addPartitionContribution( const size_t i1, const size_t j1
, const size_t i2, const size_t j2
, const Z_type partFunct
, const bool isHybridZ
, Z_type & partZ_noED );

/**
* Calls updateOptima() for each buffered entry of Z_partition. This hook is
* still called by reportOptima(), but complete tracker-free partitions were
* already consumed by updateZ() and are therefore not present in the map.
*/
virtual
void
Expand Down
2 changes: 1 addition & 1 deletion src/IntaRNA/PredictorMfeEns2d.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -194,7 +194,7 @@ fillHybridZ( const size_t j1, const size_t j2

// update mfe if needed
if (callUpdateZ) {
updateZ(i1, j1, i2, j2, curZ, true);
updateCompleteZ(i1, j1, i2, j2, curZ, true);
}

} // complementary base pair
Expand Down
5 changes: 3 additions & 2 deletions src/IntaRNA/PredictorMfeEns2d.h
Original file line number Diff line number Diff line change
Expand Up @@ -72,13 +72,14 @@ class PredictorMfeEns2d: public PredictorMfeEns {

/**
* Computes all entries of the hybridE matrix for interactions ending in
* p=j1 and q=j2 and report all valid interactions to updateZ()
* p=j1 and q=j2 and reports every complete boundary via updateCompleteZ(),
* which invokes the virtual updateZ() hook.
*
* @param j1 end of the interaction within seq 1
* @param j2 end of the interaction within seq 2
* @param i1init smallest value for i1
* @param i2init smallest value for i2
* @param callUpdateZ whether or not updateZ() is to be called
* @param callUpdateZ whether or not updateCompleteZ() is to be called
*
*/
virtual
Expand Down
158 changes: 158 additions & 0 deletions tests/PredictorMfeEnsRegression_test.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -14,6 +14,7 @@
#include "IntaRNA/SeedHandlerNoBulge.h"

#include <cmath>
#include <stdexcept>

using namespace IntaRNA;

Expand All @@ -32,6 +33,79 @@ class InspectableEnsemblePredictor : public PredictorType {
}
};

class InspectableExactEnsemblePredictor : public PredictorMfeEns2d {
public:
InspectableExactEnsemblePredictor(
const InteractionEnergy & energy,
OutputHandler & output,
PredictionTracker * tracker)
: PredictorMfeEns2d(energy, output, tracker)
{}

size_t getPartitionCount() const {
return Z_partition.size();
}
};

class CountingPredictionTracker : public PredictionTracker {
public:
explicit CountingPredictionTracker(size_t & calls)
: calls(calls)
{}

void updateOptimumCalled(const size_t, const size_t,
const size_t, const size_t, const E_type) override
{
++calls;
}

private:
size_t & calls;
};

class InterceptingExactEnsemblePredictor : public PredictorMfeEns2d {
public:
InterceptingExactEnsemblePredictor(
const InteractionEnergy & energy,
OutputHandler & output)
: PredictorMfeEns2d(energy, output, NULL)
, updateCalls(0)
, throwOnUpdate(false)
{}

size_t getPartitionCount() const {
return Z_partition.size();
}

void initializeForUpdateTest() {
initOptima();
initZ();
}

void addCompleteForUpdateTest() {
updateCompleteZ(0, 0, 0, 0, Z_type(1), true);
}

void addBufferedForUpdateTest() {
PredictorMfeEns::updateZ(0, 0, 0, 0, Z_type(1), true);
}

size_t updateCalls;
bool throwOnUpdate;

protected:
void updateZ(const size_t i1, const size_t j1,
const size_t i2, const size_t j2,
const Z_type partZ, const bool isHybridZ) override
{
++updateCalls;
if (throwOnUpdate) {
throw std::runtime_error("intercepted complete partition");
}
PredictorMfeEns::updateZ(i1, j1, i2, j2, partZ, isHybridZ);
}
};

TEST_CASE("ensemble predictor regressions", "[PredictorMfeEns]") {

#include "testEasyLoggingSetup.icc"
Expand Down Expand Up @@ -135,4 +209,88 @@ TEST_CASE("ensemble predictor regressions", "[PredictorMfeEns]") {
REQUIRE(heuristic.getPartitionCount() == 0);
REQUIRE(heuristic.getZall() == 0.0);
}

SECTION("exact complete boundaries stream unless tracker ordering is required") {
RnaSequence target("target", "GGGG");
RnaSequence query("query", "CCCC");
AccessibilityDisabled targetAcc(target, 0, NULL);
AccessibilityDisabled queryAcc(query, 0, NULL);
ReverseAccessibility reverseQueryAcc(queryAcc);
InteractionEnergyBasePair energy(targetAcc, reverseQueryAcc, 2, 2);

for (int overlapValue = OutputConstraint::OVERLAP_NONE;
overlapValue <= OutputConstraint::OVERLAP_BOTH; ++overlapValue)
{
OutputConstraint constraint(5,
static_cast<OutputConstraint::ReportOverlap>(overlapValue),
E_INF, E_INF, false, false, false, true, false);

OutputHandlerInteractionList streamedOut(constraint, 5);
InspectableExactEnsemblePredictor streamed(energy, streamedOut, NULL);
streamed.predict();
REQUIRE(streamed.getPartitionCount() == 0);

size_t trackerCalls = 0;
OutputHandlerInteractionList bufferedOut(constraint, 5);
InspectableExactEnsemblePredictor buffered(energy, bufferedOut,
new CountingPredictionTracker(trackerCalls));
buffered.predict();
REQUIRE(buffered.getPartitionCount() > 0);
REQUIRE(trackerCalls == buffered.getPartitionCount());

REQUIRE(streamed.getZall() == buffered.getZall());
REQUIRE(streamedOut.getZ() == bufferedOut.getZ());
REQUIRE(streamedOut.reported() == bufferedOut.reported());

auto streamedInteraction = streamedOut.begin();
auto bufferedInteraction = bufferedOut.begin();
for (; streamedInteraction != streamedOut.end()
&& bufferedInteraction != bufferedOut.end();
++streamedInteraction, ++bufferedInteraction)
{
REQUIRE(**streamedInteraction == **bufferedInteraction);
}
REQUIRE(streamedInteraction == streamedOut.end());
REQUIRE(bufferedInteraction == bufferedOut.end());
}
}

SECTION("exact streaming retains the virtual update hook") {
RnaSequence target("target", "GGGG");
RnaSequence query("query", "CCCC");
AccessibilityDisabled targetAcc(target, 0, NULL);
AccessibilityDisabled queryAcc(query, 0, NULL);
ReverseAccessibility reverseQueryAcc(queryAcc);
InteractionEnergyBasePair energy(targetAcc, reverseQueryAcc, 2, 2);
OutputConstraint constraint(5, OutputConstraint::OVERLAP_BOTH,
E_INF, E_INF, false, false, false, true, false);
OutputHandlerInteractionList output(constraint, 5);
InterceptingExactEnsemblePredictor predictor(energy, output);

predictor.predict();

REQUIRE(predictor.updateCalls > 0);
REQUIRE(predictor.getPartitionCount() == 0);
}

SECTION("complete update mode is restored after an override throws") {
RnaSequence target("target", "G");
RnaSequence query("query", "C");
AccessibilityDisabled targetAcc(target, 0, NULL);
AccessibilityDisabled queryAcc(query, 0, NULL);
ReverseAccessibility reverseQueryAcc(queryAcc);
InteractionEnergyBasePair energy(targetAcc, reverseQueryAcc);
OutputConstraint constraint(1, OutputConstraint::OVERLAP_BOTH,
E_INF, E_INF, false, false, false, true, false);
OutputHandlerInteractionList output(constraint, 1);
InterceptingExactEnsemblePredictor predictor(energy, output);
predictor.initializeForUpdateTest();
predictor.throwOnUpdate = true;

REQUIRE_THROWS_AS(predictor.addCompleteForUpdateTest(), std::runtime_error);

predictor.throwOnUpdate = false;
predictor.addBufferedForUpdateTest();
REQUIRE(predictor.getPartitionCount() == 1);
}
}
Loading