Skip to content

Commit

Permalink
Merge pull request #1360 from virtualcell/dan-ss-bound-transition
Browse files Browse the repository at this point in the history 
Dan ss bound transition
  • Loading branch information
@danv61
danv61 authored Sep 26, 2024
2 parents a4e4bdb + 6605638 commit 0995c9e
Show file tree
Hide file tree
Showing 5 changed files with 192 additions and 85 deletions.
130 changes: 90 additions & 40 deletions vcell-core/src/main/java/cbit/vcell/mapping/ReactionRuleSpec.java
View file
Original file line number Diff line number Diff line change
Expand Up @@ -450,11 +450,15 @@ public TransitionCondition getTransitionCondition(Map<String, Object> analysisRe
return null;
} else if(mtpReactantList.size() == 2 && mtpProductList.size() == 2) { // we look for transition with "Bound" condition, means 2 molecules
MolecularTypePattern mtpTransitionReactant = (MolecularTypePattern)analysisResults.get(MtpReactantState + "1"); // the transitioning molecule
int totalExistsBonds = 0;
int totalExplicitStates = 0;
int numExplicitStates[] = { 0, 0};
int numExistsBonds[] = { 0, 0 };
int numSpecifiedBonds[] = { 0, 0 };
MolecularComponentPattern[] mcpExplicitStates = new MolecularComponentPattern[2]; // the site with the explicit state
MolecularComponentPattern[] mcpSpecifiedBonds = new MolecularComponentPattern[2]; // the site with the specified bond

for(int i=0; i < mtpReactantList.size(); i++) {
// the transitioning site must have explicit bond in both reactant and product (which obviously must point to the condition molecule)
// there must be only one transitioning state, the one in the transitioning molecule
MolecularTypePattern mtpReactant = mtpReactantList.get(i);
List<MolecularComponentPattern> mcpList = mtpReactant.getComponentPatternList();
if(mcpList == null || mcpList.isEmpty()) {
Expand All @@ -466,48 +470,33 @@ public TransitionCondition getTransitionCondition(Map<String, Object> analysisRe
return null; // all sites must have at least one state
}
BondType bondTypeReactant = mcpReactant.getBondType();
if(BondType.Exists == bondTypeReactant) {
if(mtpTransitionReactant == mtpReactant) {
// the transition site must be in one molecule, the condition bound site must be in the other
return null;
}
numExistsBonds[i]++; // the condition binding site has bond type "Exists"
totalExistsBonds++;
if(BondType.Specified == bondTypeReactant) { // we must have exactly one Specified bond, all the others must be Possible
numSpecifiedBonds[i]++; // the condition binding site has bond type "Exists"
mcpSpecifiedBonds[i] = mcpReactant;
} else if(BondType.Possible != bondTypeReactant) {
return null; // all other bonds in the reactant must be of type "Possible"
// we know that the 2 molecules are bound somehow, because they are in the same species
// but we don't care between which sites the bond is
}
if(cspReactant.isAny()) {
continue;
} else {
// the site with explicit state must be the same that has bond type exists
if(BondType.Exists != bondTypeReactant && mtpTransitionReactant != mtpReactant) {
return null;
}

if(!cspReactant.isAny()) { // we also need to reack explicit states
numExplicitStates[i]++;
totalExplicitStates++;
mcpExplicitStates[i] = mcpReactant;
}
numExplicitStates[i]++;
totalExplicitStates++;
}
}

if(totalExistsBonds != 1) {
return null; // the condition bound site is the only one with BondType "exists"

// the explicit state if exists, must only be on the bound site
if(mcpExplicitStates[0] != null && mcpExplicitStates[0] != mcpSpecifiedBonds[0]) {
return null;
}
if(totalExplicitStates < 1 || totalExplicitStates > 2) {
// transitioning site must be in explicit state
// the condition binding site may be explicit or any
// all the other sites must be in "Any" state
if(mcpExplicitStates[1] != null && mcpExplicitStates[1] != mcpSpecifiedBonds[1]) {
return null;
}
if(numExplicitStates[0] == 1 && numExplicitStates[1] == 1) {
// above we made sure that the transition molecule and the condition molecule are not the same
return TransitionCondition.BOUND; // must have exactly one explicit state in each molecule
} else if(numExplicitStates[0] == 0 && numExplicitStates[1] == 1 && numExistsBonds[0] == 1) {
// this else-if and the next: the only explicit site is the one transitioning
// the conditional site type may also be in the "any" state but must have bond type "Exists"

// the only situations that correctly describe condition bound transition reaction
if(numSpecifiedBonds[0] == 1 && numSpecifiedBonds[1] == 1 && numExplicitStates[0] == 1 && numExplicitStates[1] == 1) {
return TransitionCondition.BOUND;
} else if(numExplicitStates[0] == 1 && numExplicitStates[1] == 0 && numExistsBonds[1] == 1) {
} else if(numSpecifiedBonds[0] == 1 && numSpecifiedBonds[1] == 1 && totalExplicitStates == 1) {
return TransitionCondition.BOUND;
}
}
Expand Down Expand Up @@ -729,8 +718,8 @@ private void writeTransitionData(StringBuilder sb, Subtype subtype, Map<String,
}
// now let's find the Site and the State of the binding condition molecule
for(MolecularComponentPattern mcpCandidate : mtpConditionReactant.getComponentPatternList()) {
if(BondType.Exists == mcpCandidate.getBondType()) {
mcpConditionReactant = mcpCandidate; // found the bond condition site, it's the one with bond type "Exists"
if(BondType.Specified == mcpCandidate.getBondType()) {
mcpConditionReactant = mcpCandidate; // found the bond condition site, it's the one with bond type "Specified"
if(!mcpConditionReactant.getComponentStatePattern().isAny()) {
stateConditionReactant = mcpConditionReactant.getComponentStatePattern().getComponentStateDefinition().getName();
}
Expand Down Expand Up @@ -891,6 +880,7 @@ public SpeciesContext getDestroyedSpecies(SpeciesContextSpec[] speciesContextSpe
private final static String SpringSaLaDMsgSubtypeMustBeIrreversible = "This reaction subtype must be irreversible. Make the reaction ireversible.";
private final static String SpringSaLaDMsgTransmembraneBinding = "Transmembrane binding not supported";
private final static String SpringSaLaDMsgSubtypeMustBeReversible = "This reaction subtype must be reversible. Make the reaction reversible and keep Kr at 0.";
private final static String SpringSaLaDMsgBoundTransitionNeedsBinding = "No Binding reaction found that would make the binding condition of this Transition reaction possible";
private static final String SpringSaLaDMsgNoComplexes = "SpringSaLaD does not accept explicit complexes in the list of reactants";

public void gatherIssues(IssueContext issueContext, List<Issue> issueList, ReactionContext rc) {
Expand Down Expand Up @@ -919,7 +909,10 @@ public void gatherIssues(IssueContext issueContext, List<Issue> issueList, React
// }
// }
if(rc.getSimulationContext().getApplicationType() == SimulationContext.Application.SPRINGSALAD) {


if(isExcluded()) {
return; // we don't care as long as it's excluded
}
List<ReactantPattern> rpList = reactionRule.getReactantPatterns();
List<ProductPattern> ppList = reactionRule.getProductPatterns();
if(rpList == null || ppList == null) {
Expand Down Expand Up @@ -950,6 +943,8 @@ public void gatherIssues(IssueContext issueContext, List<Issue> issueList, React
return;
}



// no reaction can be on the membrane
// dan 7/31/24 correction: a reaction involving a Membrane molecule can be on a membrane
// String reactionStructure = reactionRule.getStructure().getName();
Expand Down Expand Up @@ -1023,7 +1018,6 @@ public void gatherIssues(IssueContext issueContext, List<Issue> issueList, React
// }

// ------------------------------------------------------------------------------------------------------------------------

Map<String, Object> analysisResults = new LinkedHashMap<> ();
analizeReaction(analysisResults);
Subtype subtype = getSubtype(analysisResults);
Expand Down Expand Up @@ -1072,9 +1066,65 @@ public void gatherIssues(IssueContext issueContext, List<Issue> issueList, React
String tip = GenericTip;
issueList.add(new Issue(r, issueContext, IssueCategory.Identifiers, msg, tip, Issue.Severity.WARNING));
}
// these reactions cannot be reversible

// only for conditional bound transition reactions: we need a binding reaction to create the conditional bond
switch(subtype) {
case TRANSITION:
TransitionCondition tc = getTransitionCondition(analysisResults);
if(tc != TransitionCondition.BOUND) {
break; // we do this only for conditional bound transition reactions
}
// we need to find the molecule / site combination that needs to be bound
MolecularType mtTransition = ((MolecularTypePattern)analysisResults.get(MtpReactantState + "1")).getMolecularType(); // the transitioning molecule
MolecularComponent mcTransition = ((MolecularComponentPattern)analysisResults.get(McpReactantState + "1")).getMolecularComponent();
Bond bond = ((MolecularComponentPattern)analysisResults.get(McpReactantState + "1")).getBond(); // we know it's a bound transition, so it must have a bond for sure
MolecularType mtCondition = bond.molecularTypePattern.getMolecularType();
MolecularComponent mcCondition = bond.molecularComponentPattern.getMolecularComponent();

List <ReactionRule> rrList = rc.getModel().getRbmModelContainer().getReactionRuleList();
reactionRule.getReactantPattern(0).getSpeciesPattern();
// we go through all reaction rules, trying to find a binding reaction that would make this
// conditional bound transition reaction possible
boolean found = false;
for(ReactionRule rr : rrList) {
ReactionRuleSpec rrs = rc.getReactionRuleSpec(rr);
if(rrs.isExcluded()) {
continue; // even if it's the right reaction, if it's excluded it won't do the job
}
Map<String, Object> ar = new LinkedHashMap<>();
rrs.analizeReaction(ar);
Subtype st = rrs.getSubtype(ar);
if (st != Subtype.BINDING) {
continue;
}
// the 2 reactants should match the pairs mtTransition / mcTransition and respectively mtCondition / mcCondition
MolecularType mt1 = ((MolecularTypePattern) ar.get(MtpReactantBond + "1")).getMolecularType(); // '1' - reactant. left side of the bond
MolecularType mt2 = ((MolecularTypePattern) ar.get(MtpReactantBond + "2")).getMolecularType(); // '2' - reactant, right side of the bond
MolecularComponent mc1 = ((MolecularComponentPattern) ar.get(McpReactantBond + "1")).getMolecularComponent();
MolecularComponent mc2 = ((MolecularComponentPattern) ar.get(McpReactantBond + "2")).getMolecularComponent();
if (mtTransition == mt1 && mcTransition == mc1 && mtCondition == mt2 && mcCondition == mc2) {
found = true; // found the binding reaction we needed
break; // exit the for loop
} else if (mtTransition == mt2 && mcTransition == mc2 && mtCondition == mt1 && mcCondition == mc1) {
found = true;
break;
}
}
if (found == false) {
// didn't find any binding reaction that would allow our conditional bound transition reaction to happen
String msg = SpringSaLaDMsgBoundTransitionNeedsBinding;
String tip = "Create a Binding reaction to make the binding condition of this Transition reaction possible.";
issueList.add(new Issue(r, issueContext, IssueCategory.Identifiers, msg, tip, Issue.Severity.WARNING));
return;
}
// TODO: to be very strict, we should also check if we have a seed species with the right site / state combination
// needed for the transition reactant, or another transition (or creation?) reaction that would produce it
break;
default:
break;
}

switch(subtype) { // these reactions cannot be reversible
case CREATION:
case DECAY:
case TRANSITION:
Expand Down
43 changes: 35 additions & 8 deletions vcell-core/src/main/java/org/vcell/model/ssld/SsldUtils.java
View file
Original file line number Diff line number Diff line change
Expand Up @@ -926,6 +926,12 @@ private void adjustTransitionReactionSites(TransitionReaction ssldReaction, Mapp
Molecule ssldConditionalMolecule = ssldReaction.getConditionalMolecule();
Molecule ssldTransitionMolecule = ssldReaction.getMolecule();

// these will be null unless we have a conditional bound transition reaction
MolecularTypePattern mtpConditionReactant;
MolecularComponentPattern mcpConditionReactant;
MolecularTypePattern mtpConditionProduct;
MolecularComponentPattern mcpConditionProduct;

int transitionMoleculeIndexInSpeciesPattern;
if(ssldConditionalMolecule != null) {
if(!ssldReaction.getCondition().equals(TransitionReaction.BOUND_CONDITION)) {
Expand All @@ -940,30 +946,35 @@ private void adjustTransitionReactionSites(TransitionReaction ssldReaction, Mapp
State ssldConditionalState = ssldReaction.getConditionalState();

SpeciesPattern spReactant = m.getConditionReactantPattern(ssldConditionalMolecule).getSpeciesPattern();
MolecularTypePattern mtpConditionReactant = spReactant.getMolecularTypePatterns(ssldConditionalMolecule.getName()).get(0);
MolecularComponentPattern mcpConditionReactant = mtpConditionReactant.getMolecularComponentPattern(conditionalType.getName());
mtpConditionReactant = spReactant.getMolecularTypePatterns(ssldConditionalMolecule.getName()).get(0);
mcpConditionReactant = mtpConditionReactant.getMolecularComponentPattern(conditionalType.getName());
if(ssldConditionalState != TransitionReaction.ANY_STATE) {
String ssldConditionalStateName = ssldConditionalState.getName();
ComponentStateDefinition csdConditionReactant = mcpConditionReactant.getMolecularComponent().getComponentStateDefinition(ssldConditionalStateName);
ComponentStatePattern cspConditionalReactant = new ComponentStatePattern(csdConditionReactant);
mcpConditionReactant.setComponentStatePattern(cspConditionalReactant);
}
mcpConditionReactant.setBondType(MolecularComponentPattern.BondType.Exists);
mcpConditionReactant.setBond(null);
// we'll build the Bond right at the end, we don't know yet
// the other end of the bond (the molecule and site of the transition reactant)

SpeciesPattern spProduct = m.getConditionProductPattern(ssldConditionalMolecule).getSpeciesPattern();
MolecularTypePattern mtpConditionProduct = spProduct.getMolecularTypePatterns(ssldConditionalMolecule.getName()).get(0);
MolecularComponentPattern mcpConditionProduct = mtpConditionProduct.getMolecularComponentPattern(conditionalType.getName());
mtpConditionProduct = spProduct.getMolecularTypePatterns(ssldConditionalMolecule.getName()).get(0);
mcpConditionProduct = mtpConditionProduct.getMolecularComponentPattern(conditionalType.getName());
if(ssldConditionalState != TransitionReaction.ANY_STATE) {
String ssldConditionalStateName = ssldConditionalState.getName();
ComponentStateDefinition csdConditionProduct = mcpConditionProduct.getMolecularComponent().getComponentStateDefinition(ssldConditionalStateName);
ComponentStatePattern cspConditionalProduct = new ComponentStatePattern(csdConditionProduct);
mcpConditionProduct.setComponentStatePattern(cspConditionalProduct);
}
mcpConditionProduct.setBondType(MolecularComponentPattern.BondType.Exists);
mcpConditionProduct.setBond(null);
// we'll build the Bond right at the end, we don't know yet
// the other end of the bond (the molecule and site of the transition product)
} else {
transitionMoleculeIndexInSpeciesPattern = 0;

mtpConditionReactant = null;
mcpConditionReactant = null;
mtpConditionProduct = null;
mcpConditionProduct = null;
}
SiteType ssldSiteType = ssldReaction.getType();
State ssldInitialState = ssldReaction.getInitialState();
Expand All @@ -984,6 +995,14 @@ private void adjustTransitionReactionSites(TransitionReaction ssldReaction, Mapp
}
mcp.setBondType(MolecularComponentPattern.BondType.None);
}
} else if(ssldReaction.getCondition().equals(TransitionReaction.BOUND_CONDITION)) {
// we know for sure that this is the only explicit bond in each of the 2 molecular type patterns of the reactant
int bondId = 1; // correct would be: sp.nextBondId();
mcpConditionReactant.setBondId(bondId); // this also sets the BondType to Specified
mcpTransitionReactant.setBondId(bondId);
// very abusive way of bypassing the fact that the Bond constructor is private
mcpConditionReactant.setBond(SpeciesPattern.generateBond(mtpTransitionReactant, mcpTransitionReactant));
mcpTransitionReactant.setBond(SpeciesPattern.generateBond(mtpConditionReactant, mcpConditionReactant));
}

SpeciesPattern spProduct = m.getProductPatternOneFromMolecule(ssldTransitionMolecule).getSpeciesPattern();
Expand All @@ -1001,6 +1020,14 @@ private void adjustTransitionReactionSites(TransitionReaction ssldReaction, Mapp
}
mcp.setBondType(MolecularComponentPattern.BondType.None);
}
} else if(ssldReaction.getCondition().equals(TransitionReaction.BOUND_CONDITION)) {
// we know for sure that this is the only explicit bond in each of the 2 molecular type patterns of the product
int bondId = 1; // correct would be: sp.nextBondId();
mcpConditionProduct.setBondId(bondId); // this also sets the BondType to Specified
mcpTransitionProduct.setBondId(bondId);
// very abusive way of bypassing the fact that the Bond constructor is private
mcpConditionProduct.setBond(SpeciesPattern.generateBond(mtpTransitionProduct, mcpTransitionProduct));
mcpTransitionProduct.setBond(SpeciesPattern.generateBond(mtpConditionProduct, mcpConditionProduct));
}
}

Expand Down
7 changes: 5 additions & 2 deletions vcell-core/src/main/java/org/vcell/solver/langevin/LangevinLngvWriter.java
View file
Original file line number Diff line number Diff line change
Expand Up @@ -660,7 +660,7 @@ private static void writeTransitionReactions(StringBuilder sb) {
// identify the binding condition site and state
for(ParticleMolecularComponentPattern pmcp : pmtpConditionReactant.getMolecularComponentPatternList()) {
// condition state may be "any". We need to look for an "existing" bond
if(!(pmcp.getBondType() == ParticleMolecularComponentPattern.ParticleBondType.Exists)) {
if(!(pmcp.getBondType() == ParticleMolecularComponentPattern.ParticleBondType.Specified)) {
continue;
}
ParticleComponentStatePattern pcsp = pmcp.getComponentStatePattern();
Expand Down Expand Up @@ -697,6 +697,9 @@ private static void writeTransitionReactions(StringBuilder sb) {
break; // found the one and only transition
}

if(pmcpConditionReactant == null) {
System.out.println("NPE on pmcpConditionReactant");
}
// finally build the transition block
sb.append("'").append(lpjp.getName()).append("' :: ");
sb.append("'").append(pmtpTransitionReactant.getMolecularType().getName()).append("' : '")
Expand All @@ -709,7 +712,7 @@ private static void writeTransitionReactions(StringBuilder sb) {
if(TransitionCondition.BOUND == transitionCondition) {
sb.append(" '").append(pmtpConditionReactant.getMolecularType().getName()).append("' : '")
.append(pmcpConditionReactant.getMolecularComponent().getName()).append("' : '")
.append(pmcpConditionReactant.getComponentStatePattern().isAny() ? ReactionRuleSpec.ANY_STATE_STRING : pcsdConditionReactant);
.append(pmcpConditionReactant.getComponentStatePattern().isAny() ? ReactionRuleSpec.ANY_STATE_STRING : pcsdConditionReactant.getName());
}
sb.append("\n");
} // end if Subtype.TRANSITION
Expand Down
Loading

0 comments on commit 0995c9e

Please sign in to comment.