diff --git a/README.org b/README.org index 5be2f709..03efb4cd 100644 --- a/README.org +++ b/README.org @@ -1,5 +1,9 @@ #+TITLE: Notes for Rasi's modification to the program +[2018-11-23 Fri] +- For my polymer neighborhood search, I was previously looking at only the molecule components that got changed by a transformation. This results in searching of only polymer neighborhoods of the changed component, but not other unchanged components which might nevertheless be constrained by the reaction rule. +- So I now look at all molecule components that are constrained as well. This is done by looking at the template molecule corresponding to the transformation reactant. + [2018-08-23 Thu] - NFsim is working robustly now with or without polymer and network connectivity options. - As far as I can tell, both these modules are working as expected. diff --git a/src/NFcore/NFcore.hh b/src/NFcore/NFcore.hh index 833c5b4c..494c9356 100644 --- a/src/NFcore/NFcore.hh +++ b/src/NFcore/NFcore.hh @@ -960,10 +960,10 @@ namespace NFcore void depthFirstSearch(vector &members); /* functions for searching a polymer molecule within the interaction distance - * of the site where they will undergo a change + * of the site that is specified or changed * Arvind Rasi Subramaniam */ - void traversePolymerNeighborhood(vector &members, Mapping * mapping); + void traversePolymerNeighborhood(vector &members, int cIndex); /* This function is essentially same as breadthFirstSearch but adapted for * product retrieval of non-polymer molecules. * Arvind Rasi Subramaniam diff --git a/src/NFcore/molecule.cpp b/src/NFcore/molecule.cpp index e8ffc919..c3126239 100644 --- a/src/NFcore/molecule.cpp +++ b/src/NFcore/molecule.cpp @@ -416,21 +416,22 @@ void Molecule::printBondDetails(ostream &o) { */ void Molecule::printBondDetails(NFstream &o) { - int degree = 0; o<< parentMoleculeType->getName() << "\t"<getComponentName(c); - if (parentMoleculeType->getComponentStateName(c,component[c]) != "NO_STATE") { - o<< "-" << parentMoleculeType->getComponentStateName(c,component[c]); + if (this->getMoleculeType()->getSystem()->getTrackConnected()) { + o<<"\t"; + for(int c=0; cgetComponentName(c); + if (parentMoleculeType->getComponentStateName(c,component[c]) != "NO_STATE") { + o<< "-" << parentMoleculeType->getComponentStateName(c,component[c]); + } + o<<":"; + o<getMoleculeType()->getComponentName(this->indexOfBond[c]); + o<<"-"<getMoleculeTypeName()<<"_"<getUniqueID(); } - o<<":"; - o<getMoleculeType()->getComponentName(this->indexOfBond[c]); - o<<"-"<getMoleculeTypeName()<<"_"<getUniqueID(); } } o.flush(); @@ -804,9 +805,8 @@ void Molecule::printMoleculeList(vector &members) * @param mapping - contains the molecule component that changed. * @author Arvind Rasi Subramaniam */ -void Molecule::traversePolymerNeighborhood(vector &members, Mapping * mapping) { +void Molecule::traversePolymerNeighborhood(vector &members, int cIndex) { - int cIndex; int nearbycIndex; int polymerType; int polymerLocation; @@ -814,7 +814,6 @@ void Molecule::traversePolymerNeighborhood(vector &members, Mapping Molecule * mol = this; MoleculeType * mt = this->getMoleculeType(); - cIndex = mapping->getIndex(); // The molecule has no bonds to check, for eg. dna in rasi's translation model if (mol->numOfComponents == 0 | cIndex == -1) return; // if molecule is not a polymer, check the binding partner diff --git a/src/NFcore/templateMolecule.cpp b/src/NFcore/templateMolecule.cpp index e7eb2365..05878778 100644 --- a/src/NFcore/templateMolecule.cpp +++ b/src/NFcore/templateMolecule.cpp @@ -451,6 +451,10 @@ void TemplateMolecule::addBond(string thisBsiteName, //If we called this, then we are adding a bond to a nonsymmetric site int compIndex=moleculeType->getCompIndexFromName(thisBsiteName); + // Add only polymeric components for getting bonded neighbors during + // polymeric search + // Arvind Rasi Subramaniam Nov 24, 2018 + if (moleculeType->getPolymerType(compIndex) > 0) specifiedComps.push_back(compIndex); //Insert the new information bondComp.push_back(compIndex); bondCompName.push_back(thisBsiteName); diff --git a/src/NFcore/templateMolecule.hh b/src/NFcore/templateMolecule.hh index bd2b2d42..cc984253 100644 --- a/src/NFcore/templateMolecule.hh +++ b/src/NFcore/templateMolecule.hh @@ -140,6 +140,10 @@ namespace NFcore bool isMoleculeTypeAndComponentPresent(MoleculeType * mt, int cIndex); + // To get all components that are specified + // Includes Bound, Empty, Constrained, Excluded and Bonded Comps + vector getAllSpecifiedComps() {return specifiedComps;}; + protected: static int TotalTemplateMoleculeCount; @@ -157,6 +161,12 @@ namespace NFcore //// be handled separately... //// 1) unique components //// 2) symmetric components + /// + + // List of all components that are specified in the TemplateMolecule + // Includes Bound, Empty, Constrained, Excluded and Bonded Comps + // Arvind Rasi Subramaniam Nov 23, 2018 + vector specifiedComps; // Which of the unique components must be empty (no bonds) diff --git a/src/NFreactions/transformations/transformationSet.cpp b/src/NFreactions/transformations/transformationSet.cpp index 9ea23400..abba5ab5 100644 --- a/src/NFreactions/transformations/transformationSet.cpp +++ b/src/NFreactions/transformations/transformationSet.cpp @@ -745,7 +745,26 @@ bool TransformationSet::getListOfProducts(MappingSet **mappingSets, // If not, you can accidentally miss searching within the interaction distance // of all site changes. // Arvind Rasi Subramaniam - molecule->traversePolymerNeighborhood(products, mappingSets[r]->get(i)); + molecule->traversePolymerNeighborhood(products, mappingSets[r]->get(i)->getIndex()); + + // Iterate over all constrained components of the corresponding TemplateMolecule + // even if these components are not changed by the transformation. + // This helps in the case of rxns such as endocleave_3_hit to update + // ribosomes that are bound but 3' to the cut site. + // In that case, the cut site and the ribosome binding site are on + // different polymers, and the ribosome binding site does not change, + // but it is constrained. + // Arvind Rasi Subramaniam Nov 23, 2018 + TemplateMolecule * tm = transformations[r].at(i)->getTemplateMolecule(); + // If there are no templatemolecules associated with the transformation + if (!tm) continue; + // Need to check neighborhood of only polymer molecules + if (!tm->getMoleculeType()->checkIfPolymer()) continue; + for (int cIndex : tm->getAllSpecifiedComps()) { + molecule->traversePolymerNeighborhood(products, cIndex); + } + + } } }