Skip to content

Commit

Permalink
Merge pull request #155 from OpenBioSim/feature_analysis
Browse files Browse the repository at this point in the history 
Add unified FEP analysis using alchemlyb
  • Loading branch information
@lohedges
lohedges authored Oct 12, 2023
2 parents 527b5f6 + 7dd1dae commit 60feb0d
Show file tree
Hide file tree
Showing 16 changed files with 1,626 additions and 326 deletions.
76 changes: 49 additions & 27 deletions nodes/playground/prepareFEP.ipynb
View file
Original file line number Diff line number Diff line change
Expand Up @@ -159,6 +159,13 @@
" BSS.Gateway.String(\n",
" help=\"The root name for the files describing the perturbation input1->input2.\"\n",
" ),\n",
")\n",
"node.addInput(\n",
" \"somd2\",\n",
" BSS.Gateway.Boolean(\n",
" help=\"Whether to generate input for SOMD2.\",\n",
" default=False,\n",
" ),\n",
")"
]
},
Expand Down Expand Up @@ -321,18 +328,28 @@
"source": [
"# Log the mapping used\n",
"writeLog(lig1, lig2, mapping)\n",
"# Are we saving output for SOMD2?\n",
"is_somd2 = node.getInput(\"somd2\")\n",
"# File root for all output.\n",
"root = node.getInput(\"output\")\n",
"BSS.IO.saveMolecules(\n",
" \"merged_at_lam0.pdb\",\n",
" merged,\n",
" \"PDB\",\n",
" {\"coordinates\": \"coordinates0\", \"bond\": \"bond0\", \"element\": \"element0\"},\n",
")\n",
"# Generate package specific input\n",
"protocol = BSS.Protocol.FreeEnergy(runtime=2 * BSS.Units.Time.femtosecond, num_lam=3)\n",
"process = BSS.Process.Somd(system1, protocol)\n",
"process.getOutput()\n",
"with zipfile.ZipFile(\"somd_output.zip\", \"r\") as zip_hnd:\n",
" zip_hnd.extractall(\".\")"
"if is_somd2:\n",
" BSS.Stream.save(system1, root)\n",
" stream_file = \"%s.bss\" % root\n",
"else:\n",
" # Generate package specific input\n",
" protocol = BSS.Protocol.FreeEnergy(\n",
" runtime=2 * BSS.Units.Time.femtosecond, num_lam=3\n",
" )\n",
" process = BSS.Process.Somd(system1, protocol)\n",
" process.getOutput()\n",
" with zipfile.ZipFile(\"somd_output.zip\", \"r\") as zip_hnd:\n",
" zip_hnd.extractall(\".\")"
]
},
{
Expand All @@ -341,12 +358,12 @@
"metadata": {},
"outputs": [],
"source": [
"root = node.getInput(\"output\")\n",
"mergedpdb = \"%s.mergeat0.pdb\" % root\n",
"pert = \"%s.pert\" % root\n",
"prm7 = \"%s.prm7\" % root\n",
"rst7 = \"%s.rst7\" % root\n",
"mapping_str = \"%s.mapping\" % root"
"if not is_somd2:\n",
" mergedpdb = \"%s.mergeat0.pdb\" % root\n",
" pert = \"%s.pert\" % root\n",
" prm7 = \"%s.prm7\" % root\n",
" rst7 = \"%s.rst7\" % root\n",
" mapping_str = \"%s.mapping\" % root"
]
},
{
Expand All @@ -355,18 +372,19 @@
"metadata": {},
"outputs": [],
"source": [
"os.replace(\"merged_at_lam0.pdb\", mergedpdb)\n",
"os.replace(\"somd.pert\", pert)\n",
"os.replace(\"somd.prm7\", prm7)\n",
"os.replace(\"somd.rst7\", rst7)\n",
"os.replace(\"somd.mapping\", mapping_str)\n",
"try:\n",
" os.remove(\"somd_output.zip\")\n",
" os.remove(\"somd.cfg\")\n",
" os.remove(\"somd.err\")\n",
" os.remove(\"somd.out\")\n",
"except Exception:\n",
" pass"
"if not is_somd2:\n",
" os.replace(\"merged_at_lam0.pdb\", mergedpdb)\n",
" os.replace(\"somd.pert\", pert)\n",
" os.replace(\"somd.prm7\", prm7)\n",
" os.replace(\"somd.rst7\", rst7)\n",
" os.replace(\"somd.mapping\", mapping_str)\n",
" try:\n",
" os.remove(\"somd_output.zip\")\n",
" os.remove(\"somd.cfg\")\n",
" os.remove(\"somd.err\")\n",
" os.remove(\"somd.out\")\n",
" except Exception:\n",
" pass"
]
},
{
Expand All @@ -375,7 +393,11 @@
"metadata": {},
"outputs": [],
"source": [
"node.setOutput(\"nodeoutput\", [mergedpdb, pert, prm7, rst7, mapping_str])"
"if is_somd2:\n",
" output = [stream_file]\n",
"else:\n",
" output = [mergedpdb, pert, prm7, rst7, mapping_str]\n",
"node.setOutput(\"nodeoutput\", output)"
]
},
{
Expand All @@ -390,7 +412,7 @@
],
"metadata": {
"kernelspec": {
"display_name": "Python 3",
"display_name": "Python 3 (ipykernel)",
"language": "python",
"name": "python3"
},
Expand All @@ -404,7 +426,7 @@
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.7.3"
"version": "3.10.11"
}
},
"nbformat": 4,
Expand Down
74 changes: 48 additions & 26 deletions nodes/playground/prepareFEP.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -49,7 +49,7 @@ def writeLog(ligA, ligB, mapping):


def loadMapping(mapping_file):
"""Parse a text file that specifies mappings between atomic indices in input1 --> atoms in input2."""
"""Parse a text file that specifies mappings between atomic indices in input1 --> atoms in input2"""
stream = open(mapping_file, "r")
buffer = stream.readlines()
stream.close()
Expand Down Expand Up @@ -133,6 +133,13 @@ def loadMapping(mapping_file):
help="The root name for the files describing the perturbation input1->input2."
),
)
node.addInput(
"somd2",
BSS.Gateway.Boolean(
help="Whether to generate input for SOMD2.",
default=False,
),
)


# In[ ]:
Expand Down Expand Up @@ -258,52 +265,67 @@ def loadMapping(mapping_file):

# Log the mapping used
writeLog(lig1, lig2, mapping)
# Are we saving output for SOMD2?
is_somd2 = node.getInput("somd2")
# File root for all output.
root = node.getInput("output")
BSS.IO.saveMolecules(
"merged_at_lam0.pdb",
merged,
"PDB",
{"coordinates": "coordinates0", "bond": "bond0", "element": "element0"},
)
# Generate package specific input
protocol = BSS.Protocol.FreeEnergy(runtime=2 * BSS.Units.Time.femtosecond, num_lam=3)
process = BSS.Process.Somd(system1, protocol)
process.getOutput()
with zipfile.ZipFile("somd_output.zip", "r") as zip_hnd:
zip_hnd.extractall(".")
if is_somd2:
BSS.Stream.save(system1, root)
stream_file = "%s.bss" % root
else:
# Generate package specific input
protocol = BSS.Protocol.FreeEnergy(
runtime=2 * BSS.Units.Time.femtosecond, num_lam=3
)
process = BSS.Process.Somd(system1, protocol)
process.getOutput()
with zipfile.ZipFile("somd_output.zip", "r") as zip_hnd:
zip_hnd.extractall(".")


# In[ ]:


root = node.getInput("output")
mergedpdb = "%s.mergeat0.pdb" % root
pert = "%s.pert" % root
prm7 = "%s.prm7" % root
rst7 = "%s.rst7" % root
mapping_str = "%s.mapping" % root
if not is_somd2:
mergedpdb = "%s.mergeat0.pdb" % root
pert = "%s.pert" % root
prm7 = "%s.prm7" % root
rst7 = "%s.rst7" % root
mapping_str = "%s.mapping" % root


# In[ ]:


os.replace("merged_at_lam0.pdb", mergedpdb)
os.replace("somd.pert", pert)
os.replace("somd.prm7", prm7)
os.replace("somd.rst7", rst7)
os.replace("somd.mapping", mapping_str)
try:
os.remove("somd_output.zip")
os.remove("somd.cfg")
os.remove("somd.err")
os.remove("somd.out")
except Exception:
pass
if not is_somd2:
os.replace("merged_at_lam0.pdb", mergedpdb)
os.replace("somd.pert", pert)
os.replace("somd.prm7", prm7)
os.replace("somd.rst7", rst7)
os.replace("somd.mapping", mapping_str)
try:
os.remove("somd_output.zip")
os.remove("somd.cfg")
os.remove("somd.err")
os.remove("somd.out")
except Exception:
pass


# In[ ]:


node.setOutput("nodeoutput", [mergedpdb, pert, prm7, rst7, mapping_str])
if is_somd2:
output = [stream_file]
else:
output = [mergedpdb, pert, prm7, rst7, mapping_str]
node.setOutput("nodeoutput", output)


# In[ ]:
Expand Down
12 changes: 9 additions & 3 deletions python/BioSimSpace/Align/_align.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -36,8 +36,6 @@
import csv as _csv
import os as _os
import subprocess as _subprocess
import shutil as _shutil
import shlex as _shlex
import sys as _sys

from .._Utils import _try_import, _have_imported, _assert_imported
Expand Down Expand Up @@ -925,7 +923,7 @@ def matchAtoms(
mcs_smarts = _Chem.MolFromSmarts(mcs.smartsString)

except:
raise RuntimeError("RDKIT MCS mapping failed!")
raise RuntimeError("RDKit MCS mapping failed!")

# Score the mappings and return them in sorted order (best to worst).
mappings, scores = _score_rdkit_mappings(
Expand All @@ -948,6 +946,14 @@ def matchAtoms(
if prematch != {}:
_warnings.warn("RDKit mapping didn't include prematch. Using Sire MCS.")

# Sire MCS isn't currently supported on Windows.
if _sys.platform == "win32":
_warnings.warn("Sire MCS is currently unsupported on Windows.")
if return_scores:
return {}, []
else:
return {}

# Warn about unsupported options.
if not complete_rings_only:
_warnings.warn(
Expand Down
Loading

0 comments on commit 60feb0d

Please sign in to comment.