1 #!/bin/env python 2 # 3 # File: OpenMMPrepareMacromolecule.py 4 # Author: Manish Sud <msud@san.rr.com> 5 # 6 # Copyright (C) 2025 Manish Sud. All rights reserved. 7 # 8 # The functionality available in this script is implemented using OpenMM, an 9 # open source molecuar simulation package. 10 # 11 # This file is part of MayaChemTools. 12 # 13 # MayaChemTools is free software; you can redistribute it and/or modify it under 14 # the terms of the GNU Lesser General Public License as published by the Free 15 # Software Foundation; either version 3 of the License, or (at your option) any 16 # later version. 17 # 18 # MayaChemTools is distributed in the hope that it will be useful, but without 19 # any warranty; without even the implied warranty of merchantability of fitness 20 # for a particular purpose. See the GNU Lesser General Public License for more 21 # details. 22 # 23 # You should have received a copy of the GNU Lesser General Public License 24 # along with MayaChemTools; if not, see <http://www.gnu.org/licenses/> or 25 # write to the Free Software Foundation Inc., 59 Temple Place, Suite 330, 26 # Boston, MA, 02111-1307, USA. 27 # 28 29 from __future__ import print_function 30 31 # Add local python path to the global path and import standard library modules... 32 import os 33 import sys; sys.path.insert(0, os.path.join(os.path.dirname(sys.argv[0]), "..", "lib", "Python")) 34 import time 35 import re 36 37 # OpenMM imports... 38 try: 39 import openmm as mm 40 import pdbfixer 41 from pdbfixer.pdbfixer import __version__ as pdbfixerversion 42 except ImportError as ErrMsg: 43 sys.stderr.write("\nFailed to import OpenMM related module/package: %s\n" % ErrMsg) 44 sys.stderr.write("Check/update your OpenMM environment and try again.\n\n") 45 sys.exit(1) 46 47 # MayaChemTools imports... 48 try: 49 from docopt import docopt 50 import MiscUtil 51 import OpenMMUtil 52 except ImportError as ErrMsg: 53 sys.stderr.write("\nFailed to import MayaChemTools module/package: %s\n" % ErrMsg) 54 sys.stderr.write("Check/update your MayaChemTools environment and try again.\n\n") 55 sys.exit(1) 56 57 ScriptName = os.path.basename(sys.argv[0]) 58 Options = {} 59 OptionsInfo = {} 60 61 def main(): 62 """Start execution of the script.""" 63 64 MiscUtil.PrintInfo("\n%s (OpenMM v%s; PDBFixerVersion v%s; MayaChemTools v%s; %s): Starting...\n" % (ScriptName, mm.Platform.getOpenMMVersion(), pdbfixerversion, MiscUtil.GetMayaChemToolsVersion(), time.asctime())) 65 66 (WallClockTime, ProcessorTime) = MiscUtil.GetWallClockAndProcessorTime() 67 68 # Retrieve command line arguments and options... 69 RetrieveOptions() 70 71 # Process and validate command line arguments and options... 72 ProcessOptions() 73 74 # Perform actions required by the script... 75 PrepareMacromolecule() 76 77 MiscUtil.PrintInfo("\n%s: Done...\n" % ScriptName) 78 MiscUtil.PrintInfo("Total time: %s" % MiscUtil.GetFormattedElapsedTime(WallClockTime, ProcessorTime)) 79 80 def PrepareMacromolecule(): 81 """Prepare a macromolecule for simulation and write it out.""" 82 83 # Read macromolecule... 84 PDBFixerObjectHandle = ReadMacromolecule() 85 86 # Identify missing residues... 87 IdentifyMissingResidues(PDBFixerObjectHandle) 88 89 # Identify and replace non-standard residues... 90 IdentifyAndReplaceNonStandardResidues(PDBFixerObjectHandle) 91 92 # Identify and delete heterogen residues... 93 IdentifyAndDeleteHeterogenResidues(PDBFixerObjectHandle) 94 95 # Identify and add missing atoms... 96 IdentifyAndAddMissingAtoms(PDBFixerObjectHandle) 97 98 if OptionsInfo["WaterBox"]: 99 AddWaterBox(PDBFixerObjectHandle) 100 elif OptionsInfo["Membrane"]: 101 AddLipidMembrane(PDBFixerObjectHandle) 102 else: 103 MiscUtil.PrintInfo("\nSkipping addtion of water box or membrane...") 104 105 # Write macromolecule... 106 WriteMacromolecule(PDBFixerObjectHandle) 107 108 def ReadMacromolecule(): 109 """Read macromolecule.""" 110 111 Infile = OptionsInfo["Infile"] 112 MiscUtil.PrintInfo("\nProcessing file %s..." % Infile) 113 114 PDBFixerObjectHandle = pdbfixer.pdbfixer.PDBFixer(filename = Infile) 115 116 return PDBFixerObjectHandle 117 118 def WriteMacromolecule(PDBFixerObjectHandle): 119 """Write macromolecule. """ 120 121 MiscUtil.PrintInfo("\nGenerating output file %s..." % OptionsInfo["Outfile"]) 122 OpenMMUtil.WritePDBFile(OptionsInfo["Outfile"], PDBFixerObjectHandle.topology, PDBFixerObjectHandle.positions, KeepIDs = True) 123 124 def IdentifyMissingResidues(PDBFixerObjectHandle): 125 """Identify missing residues based on PDB REMARK records.""" 126 127 MiscUtil.PrintInfo("\nIdentifying missing residues...") 128 129 PDBFixerObjectHandle.missingResidues = {} 130 PDBFixerObjectHandle.findMissingResidues() 131 132 ListMissingResidues(PDBFixerObjectHandle) 133 134 def ListMissingResidues(PDBFixerObjectHandle): 135 """List missing residues.""" 136 137 if len(PDBFixerObjectHandle.missingResidues) == 0: 138 MiscUtil.PrintInfo("Number of missing residues: 0") 139 return 140 141 Chains = list(PDBFixerObjectHandle.topology.chains()) 142 143 MissingResiduesCount = 0 144 MissingResiduesInfo = [] 145 for Index, Key in enumerate(sorted(PDBFixerObjectHandle.missingResidues)): 146 ChainIndex = Key[0] 147 # ResidueIndex corresponds to the residue after the insertion of missing residues...` 148 ResidueIndex = Key[1] 149 150 MissingResidues = PDBFixerObjectHandle.missingResidues[Key] 151 MissingResiduesCount += len(MissingResidues) 152 153 Chain = Chains[ChainIndex] 154 ChainResidues = list(Chain.residues()) 155 156 if ResidueIndex < len(ChainResidues): 157 ResidueNumOffset = int(ChainResidues[ResidueIndex].id) - len(MissingResidues) - 1 158 else: 159 ResidueNumOffset = int(ChainResidues[-1].id) 160 161 StartResNum = ResidueNumOffset + 1 162 EndResNum = ResidueNumOffset + len(MissingResidues) 163 164 MissingResiduesInfo.append("Chain: %s; StartResNum-EndResNum: %s-%s; Count: %s; Residues: %s" % (Chain.id, StartResNum, EndResNum, len(MissingResidues), ", ".join(MissingResidues))) 165 166 MiscUtil.PrintInfo("Total number of missing residues: %s" % (MissingResiduesCount)) 167 if OptionsInfo["ListDetails"]: 168 MiscUtil.PrintInfo("%s" % ("\n".join(MissingResiduesInfo))) 169 170 def IdentifyAndReplaceNonStandardResidues(PDBFixerObjectHandle): 171 """Identify and replace non-standard residues.""" 172 173 MiscUtil.PrintInfo("\nIdentifying non-standard residues...") 174 175 PDBFixerObjectHandle.findNonstandardResidues() 176 NonStandardResidues = PDBFixerObjectHandle.nonstandardResidues 177 NumNonStandardResiues = len(NonStandardResidues) 178 179 MiscUtil.PrintInfo("Total number of non-standard residues: %s" % NumNonStandardResiues) 180 181 if NumNonStandardResiues == 0: 182 return 183 184 if OptionsInfo["ListDetails"]: 185 ListNonStandardResidues(PDBFixerObjectHandle) 186 187 if not OptionsInfo["ReplaceNonStandardResidues"]: 188 MiscUtil.PrintInfo("Skipping replacement of non-standard residues...") 189 return 190 191 MiscUtil.PrintInfo("Replacing non-standard residues...") 192 PDBFixerObjectHandle.replaceNonstandardResidues() 193 194 def ListNonStandardResidues(PDBFixerObjectHandle): 195 """List non-standard residues. """ 196 197 NonStandardResidues = PDBFixerObjectHandle.nonstandardResidues 198 199 if len(NonStandardResidues) == 0: 200 return 201 202 NonStandardResiduesMappingInfo = [] 203 for Values in NonStandardResidues: 204 NonStandardRes = Values[0] 205 StandardResName = Values[1] 206 207 MappingInfo = "<%s %s %s>: %s" % (NonStandardRes.id, NonStandardRes.name, NonStandardRes.chain.id, StandardResName) 208 NonStandardResiduesMappingInfo.append(MappingInfo) 209 210 if len(NonStandardResiduesMappingInfo): 211 MiscUtil.PrintInfo("Non-standard residues mapping: %s\n" % ",".join(NonStandardResiduesMappingInfo)) 212 213 def IdentifyAndDeleteHeterogenResidues(PDBFixerObjectHandle): 214 """Identify and delete heterogen residues.""" 215 216 DeleteHeterogens = OptionsInfo["DeleteHeterogens"] 217 if DeleteHeterogens is None: 218 MiscUtil.PrintInfo("\nSkipping deletion of any heterogen residues...") 219 return 220 221 if re.match("^All$", DeleteHeterogens, re.I): 222 MiscUtil.PrintInfo("\nDeleting all heterogen residues including water...") 223 PDBFixerObjectHandle.removeHeterogens(keepWater = False) 224 elif re.match("^AllExceptWater$", DeleteHeterogens, re.I): 225 MiscUtil.PrintInfo("\nDeleting all heterogen residues except water...") 226 PDBFixerObjectHandle.removeHeterogens(keepWater = True) 227 elif re.match("^WaterOnly$", DeleteHeterogens, re.I): 228 MiscUtil.PrintInfo("\nDeleting water only heterogen residues...") 229 DeleteWater(PDBFixerObjectHandle) 230 else: 231 MiscUtil.PrintError("The value specified, %s, for option \"--deleteHeterogens\" is not valid. Supported values: All AllExceptWater WaterOnly None" % DeleteHeterogens) 232 233 def DeleteWater(PDBFixerObjectHandle): 234 """Delete water. """ 235 236 ModellerObjectHandle = mm.app.Modeller(PDBFixerObjectHandle.topology, PDBFixerObjectHandle.positions) 237 ModellerObjectHandle.deleteWater() 238 239 PDBFixerObjectHandle.topology = ModellerObjectHandle.topology 240 PDBFixerObjectHandle.positions = ModellerObjectHandle.positions 241 242 def IdentifyAndAddMissingAtoms(PDBFixerObjectHandle): 243 """Identify and missing atoms along with already identified missing residues.""" 244 245 MiscUtil.PrintInfo("\nIdentifying missing atoms...") 246 PDBFixerObjectHandle.findMissingAtoms() 247 248 ListMissingAtoms(PDBFixerObjectHandle) 249 250 if OptionsInfo["AddHeavyAtoms"] and OptionsInfo["AddResidues"]: 251 Msg = "Adding missing heavy atoms along with missing residues..." 252 elif OptionsInfo["AddHeavyAtoms"]: 253 Msg = "Adding missing heavy atoms..." 254 elif OptionsInfo["AddResidues"]: 255 Msg = "Adding missing residues..." 256 else: 257 Msg = "Skipping addition of any heavy atoms along with any missing residues..." 258 259 if not OptionsInfo["AddHeavyAtoms"]: 260 PDBFixerObjectHandle.missingAtoms = {} 261 PDBFixerObjectHandle.missingTerminals = {} 262 263 MiscUtil.PrintInfo("\n%s" % Msg) 264 if OptionsInfo["AddHeavyAtoms"] or OptionsInfo["AddResidues"]: 265 PDBFixerObjectHandle.addMissingAtoms() 266 267 if OptionsInfo["AddHydrogens"]: 268 MiscUtil.PrintInfo("Adding missing hydrogens at pH %s..." % OptionsInfo["AddHydrogensAtpH"]) 269 PDBFixerObjectHandle.addMissingHydrogens(pH = OptionsInfo["AddHydrogensAtpH"], forcefield = None) 270 else: 271 MiscUtil.PrintInfo("Skipping addition of any missing hydrogens...") 272 273 def ListMissingAtoms(PDBFixerObjectHandle): 274 """List missing atoms.""" 275 276 if len(PDBFixerObjectHandle.missingAtoms) == 0 and len(PDBFixerObjectHandle.missingTerminals) == 0: 277 MiscUtil.PrintInfo("Total number of missing atoms: 0") 278 return 279 280 ListMissingAtomsInfo(PDBFixerObjectHandle.missingAtoms, TerminalAtoms = False) 281 ListMissingAtomsInfo(PDBFixerObjectHandle.missingTerminals, TerminalAtoms = True) 282 283 def ListMissingAtomsInfo(MissingAtoms, TerminalAtoms): 284 """Setup missing atoms information. """ 285 286 MissingAtomsInfo = [] 287 MissingAtomsCount = 0 288 MissingAtomsResiduesCount = 0 289 for Residue, Atoms in MissingAtoms.items(): 290 MissingAtomsResiduesCount += 1 291 MissingAtomsCount += len(Atoms) 292 293 if TerminalAtoms: 294 AtomsNames = [AtomName for AtomName in Atoms] 295 else: 296 AtomsNames = [Atom.name for Atom in Atoms] 297 AtomsInfo = "<%s %s %s>: %s" % (Residue.id, Residue.name, Residue.chain.id, ", ".join(AtomsNames)) 298 MissingAtomsInfo.append(AtomsInfo) 299 300 AtomsLabel = "terminal atoms" if TerminalAtoms else "atoms" 301 if MissingAtomsCount == 0: 302 MiscUtil.PrintInfo("Total number of missing %s: %s" % (MissingAtomsCount, AtomsLabel)) 303 else: 304 ResiduesLabel = "residue" if MissingAtomsResiduesCount == 1 else "residues" 305 MiscUtil.PrintInfo("Total number of missing %s across %s %s: %s" % (AtomsLabel, MissingAtomsResiduesCount, ResiduesLabel, MissingAtomsCount)) 306 if OptionsInfo["ListDetails"]: 307 MiscUtil.PrintInfo("Missing %s across residues: %s" % (AtomsLabel, "; ".join(MissingAtomsInfo))) 308 309 def AddWaterBox(PDBFixerObjectHandle): 310 """Add water box. """ 311 312 MiscUtil.PrintInfo("\nAdding water box...") 313 314 WaterBoxParams = OptionsInfo["WaterBoxParams"] 315 316 Size, Padding, Shape = [None] * 3 317 if WaterBoxParams["ModeSize"]: 318 SizeList = WaterBoxParams["SizeList"] 319 Size = mm.Vec3(SizeList[0], SizeList[1], SizeList[2]) * mm.unit.nanometer 320 elif WaterBoxParams["ModePadding"]: 321 Padding = WaterBoxParams["Padding"] * mm.unit.nanometer 322 Shape = WaterBoxParams["Shape"] 323 else: 324 MiscUtil.PrintError("The parameter value, %s, specified for parameter name, mode, using \"--waterBoxParams\" option is not a valid value. Supported values: Size Padding \n" % (WaterBoxParams["Mode"])) 325 326 IonicStrength = OptionsInfo["IonicStrength"] * mm.unit.molar 327 PDBFixerObjectHandle.addSolvent(boxSize = Size, padding = Padding, boxShape = Shape, positiveIon = OptionsInfo["IonPositive"], negativeIon = OptionsInfo["IonNegative"], ionicStrength = IonicStrength) 328 329 def AddLipidMembrane(PDBFixerObjectHandle): 330 """Add lipid membrane along with water.""" 331 332 MiscUtil.PrintInfo("\nAdding membrane along with water...") 333 334 MembraneParams = OptionsInfo["MembraneParams"] 335 336 LipidType = MembraneParams["LipidType"] 337 MembraneCenterZ = MembraneParams["MembraneCenterZ"] * mm.unit.nanometer 338 Padding = MembraneParams["Padding"] * mm.unit.nanometer 339 340 IonicStrength = OptionsInfo["IonicStrength"] * mm.unit.molar 341 342 PDBFixerObjectHandle.addMembrane(lipidType = LipidType, membraneCenterZ = MembraneCenterZ, minimumPadding = Padding, positiveIon = OptionsInfo["IonPositive"], negativeIon = OptionsInfo["IonNegative"], ionicStrength = IonicStrength) 343 344 def ProcessMembraneParamsOption(): 345 """Process option for membrane parameters.""" 346 347 ParamsOptionName = "--membraneParams" 348 ParamsOptionValue = Options[ParamsOptionName] 349 ParamsDefaultInfo = {"LipidType" : ["str", "POPC"], "MembraneCenterZ": ["float", 0.0], "Padding": ["float", 1.0]} 350 MembraneParams = MiscUtil.ProcessOptionNameValuePairParameters(ParamsOptionName, ParamsOptionValue, ParamsDefaultInfo) 351 352 for ParamName in ["LipidType", "MembraneCenterZ", "Padding"]: 353 ParamValue = MembraneParams[ParamName] 354 if ParamName == "LipidType": 355 LipidTypes = ["POPC", "POPE", "DLPC", "DLPE", "DMPC", "DOPC", "DPPC"] 356 if ParamValue not in LipidTypes: 357 MiscUtil.PrintError("The parameter value, %s, specified for parameter name, %s, using \"%s\" option is not a valid value. Supported values: %s \n" % (ParamValue, ParamName, ParamsOptionName, LipidTypes)) 358 elif ParamName == "Padding": 359 if ParamValue <= 0: 360 MiscUtil.PrintError("The parameter value, %s, specified for parameter name, %s, using \"%s\" option is not a valid value. Supported values: > 0\n" % (ParamValue, ParamName, ParamsOptionName)) 361 362 OptionsInfo["MembraneParams"] = MembraneParams 363 364 def ProcessOptions(): 365 """Process and validate command line arguments and options.""" 366 367 MiscUtil.PrintInfo("Processing options...") 368 369 # Validate options... 370 ValidateOptions() 371 372 OptionsInfo["Infile"] = Options["--infile"] 373 FileDir, FileName, FileExt = MiscUtil.ParseFileName(OptionsInfo["Infile"]) 374 OptionsInfo["InfileRoot"] = FileName 375 376 OptionsInfo["Outfile"] = Options["--outfile"] 377 378 OptionsInfo["AddHeavyAtoms"] = True if re.match("^yes$", Options["--addHeavyAtoms"], re.I) else False 379 380 OptionsInfo["AddHydrogens"] = True if re.match("^yes$", Options["--addHydrogens"], re.I) else False 381 OptionsInfo["AddHydrogensAtpH"] = float(Options["--addHydrogensAtpH"]) 382 383 OptionsInfo["AddResidues"] = True if re.match("^yes$", Options["--addResidues"], re.I) else False 384 385 DeleteHeterogens = Options["--deleteHeterogens"] 386 if re.match("^All$", DeleteHeterogens, re.I): 387 DeleteHeterogens = "All" 388 elif re.match("^AllExceptWater$", DeleteHeterogens, re.I): 389 DeleteHeterogens = "AllExceptWater" 390 elif re.match("^WaterOnly$", DeleteHeterogens, re.I): 391 DeleteHeterogens = "WaterOnly" 392 elif re.match("^None$", DeleteHeterogens, re.I): 393 DeleteHeterogens = None 394 elif re.match("^auto$", DeleteHeterogens, re.I): 395 DeleteHeterogens = None 396 if re.match("^yes$", Options["--membrane"], re.I) or re.match("^yes$", Options["--waterBox"], re.I): 397 DeleteHeterogens = "WaterOnly" 398 else: 399 MiscUtil.PrintError("The value specified, %s, for option \"--deleteHeterogens\" is not valid. Supported values: All AllExceptWater WaterOnly None" % DeleteHeterogens) 400 OptionsInfo["DeleteHeterogens"] = DeleteHeterogens 401 402 OptionsInfo["IonPositive"] = Options["--ionPositive"] 403 OptionsInfo["IonNegative"] = Options["--ionNegative"] 404 OptionsInfo["IonicStrength"] = float(Options["--ionicStrength"]) 405 406 OptionsInfo["ListDetails"] = True if re.match("^yes$", Options["--listDetails"], re.I) else False 407 408 OptionsInfo["Membrane"] = True if re.match("^yes$", Options["--membrane"], re.I) else False 409 ProcessMembraneParamsOption() 410 411 OptionsInfo["ReplaceNonStandardResidues"] = True if re.match("^yes$", Options["--replaceNonStandardResidues"], re.I) else False 412 413 OptionsInfo["WaterBox"] = True if re.match("^yes$", Options["--waterBox"], re.I) else False 414 OptionsInfo["WaterBoxParams"] = OpenMMUtil.ProcessOptionOpenMMWaterBoxParameters("--waterBoxParams", Options["--waterBoxParams"]) 415 416 OptionsInfo["Overwrite"] = Options["--overwrite"] 417 418 def RetrieveOptions(): 419 """Retrieve command line arguments and options.""" 420 421 # Get options... 422 global Options 423 Options = docopt(_docoptUsage_) 424 425 # Set current working directory to the specified directory... 426 WorkingDir = Options["--workingdir"] 427 if WorkingDir: 428 os.chdir(WorkingDir) 429 430 # Handle examples option... 431 if "--examples" in Options and Options["--examples"]: 432 MiscUtil.PrintInfo(MiscUtil.GetExamplesTextFromDocOptText(_docoptUsage_)) 433 sys.exit(0) 434 435 def ValidateOptions(): 436 """Validate option values.""" 437 438 MiscUtil.ValidateOptionFilePath("-i, --infile", Options["--infile"]) 439 MiscUtil.ValidateOptionFileExt("-i, --infile", Options["--infile"], "pdb cif") 440 441 MiscUtil.ValidateOptionFileExt("-o, --outfile", Options["--outfile"], "pdb cif") 442 MiscUtil.ValidateOptionsOutputFileOverwrite("-o, --outfile", Options["--outfile"], "--overwrite", Options["--overwrite"]) 443 MiscUtil.ValidateOptionsDistinctFileNames("-i, --infile", Options["--infile"], "-o, --outfile", Options["--outfile"]) 444 445 MiscUtil.ValidateOptionTextValue("--addHeavyAtoms", Options["--addHeavyAtoms"], "yes no") 446 447 MiscUtil.ValidateOptionTextValue("--addHydrogens", Options["--addHydrogens"], "yes no") 448 MiscUtil.ValidateOptionFloatValue("--addHydrogensAtpH", Options["--addHydrogensAtpH"], {">=": 0}) 449 450 MiscUtil.ValidateOptionTextValue("--addResidues", Options["--addResidues"], "yes no") 451 452 MiscUtil.ValidateOptionTextValue("--deleteHeterogens", Options["--deleteHeterogens"], "All AllExceptWater WaterOnly None auto") 453 DeleteHeterogens = Options["--deleteHeterogens"] 454 if re.match("^(AllExceptWater|None)$", DeleteHeterogens, re.I): 455 if re.match("^yes$", Options["--membrane"], re.I): 456 MiscUtil.PrintError("The value specified, %s, for option \"--deleteHeterogens\" is not valid during \"yes\" value of option \"--membrane\". " % DeleteHeterogens) 457 elif re.match("^yes$", Options["--waterBox"], re.I): 458 MiscUtil.PrintError("The value specified, %s, for option \"--deleteHeterogens\" is not valid during \"yes\" value of option \"--waterBox\". " % DeleteHeterogens) 459 460 ValidValues = "Li+ Na+ K+ Rb+ Cs+" 461 EscapedValidValuesPattern = "Li\+|Na\+|K\+|Rb\+|Cs\+" 462 Value = Options["--ionPositive"] 463 if not re.match("^(%s)$" % EscapedValidValuesPattern, Value): 464 MiscUtil.PrintError("The value specified, %s, for option \"-ionPositive\" is not valid: Supported value(s): %s" % (Value, ValidValues)) 465 466 ValidValues = "F- Cl- Br- I-" 467 ValidValuesPattern = "F-|Cl-|Br-|I-" 468 Value = Options["--ionNegative"] 469 if not re.match("^(%s)$" % ValidValuesPattern, Value): 470 MiscUtil.PrintError("The value specified, %s, for option \"-ionNegative\" is not valid: Supported value(s): %s" % (Value, ValidValues)) 471 472 MiscUtil.ValidateOptionFloatValue("--ionicStrength", Options["--ionicStrength"], {">=": 0}) 473 474 MiscUtil.ValidateOptionTextValue("--listDetails", Options["--listDetails"], "yes no") 475 MiscUtil.ValidateOptionTextValue("--replaceNonStandardResidues", Options["--replaceNonStandardResidues"], "yes no") 476 477 MiscUtil.ValidateOptionTextValue("--membrane", Options["--membrane"], "yes no") 478 MiscUtil.ValidateOptionTextValue("--waterBox", Options["--waterBox"], "yes no") 479 480 if re.match("^yes$", Options["--membrane"], re.I) and re.match("^yes$", Options["--waterBox"], re.I): 481 MiscUtil.PrintError("You've specified, \"Yes\" for both \"--membrane\" and \"--waterBox\" options. It's not allowed. You must choose between adding a water box or a membrane. The water is automatically added during the construction of the membrane.") 482 483 # Setup a usage string for docopt... 484 _docoptUsage_ = """ 485 OpenMMPrepareMacromolecule.py - Prepare a macromolecule for simulation. 486 487 Usage: 488 OpenMMPrepareMacromolecule.py [--addHeavyAtoms <yes or no>] [--addHydrogens <yes or no>] [--addHydrogensAtpH <number>] 489 [--addResidues <yes or no>] [--deleteHeterogens <All, AllExceptWater, WaterOnly, None>] [--ionPositive <text>] 490 [--ionNegative <text>] [--ionicStrength <number>] [--listDetails <yes or no> ] [--membrane <yes or no>] 491 [--membraneParams <Name,Value,..>] [--overwrite] [--replaceNonStandardResidues <yes or no>] [--waterBox <yes or no>] 492 [--waterBoxParams <Name,Value,..>] [-w <dir>] -i <infile> -o <outfile> 493 OpenMMPrepareMacromolecule.py -h | --help | -e | --examples 494 495 Description: 496 Prepare a macromolecute in an input file for molecular simulation and 497 and write it out to an output file. The macromolecule is prepared by 498 automatically performing the following actions: 499 500 . Identify and replace non-standard residues 501 . Add missing residues 502 . Add missing heavy atoms 503 . Add missing hydrogens 504 505 You may optionally remove heterogens, add a water box, and add a lipid 506 membrane. 507 508 The supported input file format are: PDB (.pdb) and CIF (.cif) 509 510 The supported output file formats are: PDB (.pdb) and CIF (.cif) 511 512 Options: 513 --addHeavyAtoms <yes or no> [default: yes] 514 Add missing non-hydrogen atoms based on the templates. The terminal atoms 515 are also added. 516 --addHydrogens <yes or no> [default: yes] 517 Add missing hydrogens at pH value specified using '-addHydrogensAtpH' 518 option. The missing hydrogens are added based on the templates and pKa 519 calculations are performed. 520 --addHydrogensAtpH <number> [default: 7.0] 521 pH value to use for adding missing hydrogens. 522 --addResidues <yes or no> [default: yes] 523 Add missing residues unidentified based on the PDB records. 524 --deleteHeterogens <All, AllExceptWater, WaterOnly, None> [default: auto] 525 Delete heterogens corresponding to non-standard names of amino acids, dna, 526 and rna along with any ligand names. 'N' and 'UNK' also consider standard 527 residues. Default value: WaterOnly during addition of WaterBox or Membrane; 528 Otherwise, None. 529 530 The 'AllExceptWater' or 'None' values are not allowed during the addition of 531 a water box or membrane. The waters must be deleted as they are explicitly 532 added during the construction of a water box and membrane. 533 -e, --examples 534 Print examples. 535 -h, --help 536 Print this help message. 537 -i, --infile <infile> 538 Input file name. 539 --ionPositive <text> [default: Na+] 540 Type of positive ion to add during the addition of a water box or membrane. 541 Possible values: Li+, Na+, K+, Rb+, or Cs+. 542 --ionNegative <text> [default: Cl-] 543 Type of negative ion to add during the addition of a water box or membrane. 544 Possible values: Cl-, Br-, F-, or I-. 545 --ionicStrength <number> [default: 0.0] 546 Total concentration of both positive and negative ions to add excluding 547 the ions added to neutralize the system during the addition of a water box 548 or a membrane. 549 -l, --listDetails <yes or no> [default: no] 550 List details about missing and non-standard residues along with residues 551 containing missing atoms. 552 --membrane <yes or no> [default: no] 553 Add lipid membrane along with a water box. The script replies on OpenMM 554 method modeller.addMebrane() to perform the task. The membrane is added 555 in the XY plane. The existing macromolecule must be oriented and positioned 556 correctly. 557 558 A word to the wise: You may want to start with a model from the Orientations 559 of Proteins in Membranes (OPM) database at http://opm.phar.umich.edu. 560 561 The size of the membrane and water box are determined by the value of 562 'padding' parameter specified using '--membraneParams' option. All atoms 563 in macromolecule are guaranteed to be at least this far from any edge of 564 the periodic box. 565 --membraneParams <Name,Value,..> [default: auto] 566 A comma delimited list of parameter name and value pairs for adding 567 a lipid membrane and water. 568 569 The supported parameter names along with their default values are 570 are shown below: 571 572 lipidType, POPC [ Possible values: POPC, POPE, DLPC, DLPE, DMPC, 573 DOPC or DPPC ] 574 membraneCenterZ, 0.0 575 padding, 1.0 576 577 A brief description of parameters is provided below: 578 579 lipidType: Type of lipid to use for constructing the membrane. 580 membraneCenterZ: Position along the Z axis of the center of 581 the membrane in nanomertes. 582 padding: Minimum padding distance to use in nanomertes. It's used 583 to determine the size of the membrane and water box. All atoms 584 in the macromolecule are guaranteed to be at least this far from 585 any edge of the periodic box. 586 587 -o, --outfile <outfile> 588 Output file name. 589 --overwrite 590 Overwrite existing files. 591 --replaceNonStandardResidues <yes or no> [default: yes] 592 Replace non-standard residue names by standard residue names based on the 593 list of non-standard residues available in pdbfixer. 594 --waterBox <yes or no> [default: no] 595 Add water box. 596 --waterBoxParams <Name,Value,..> [default: auto] 597 A comma delimited list of parameter name and value pairs for adding 598 a water box. 599 600 The supported parameter names along with their default values are 601 are shown below: 602 603 mode, Padding [ Possible values: Size or Padding ] 604 size, None [ Possible value: xsize ysize zsize ] 605 padding, 1.0 606 shape, cube [ Possible values: cube, dodecahedron, or octahedron ] 607 608 A brief description of parameters is provided below: 609 610 mode: Specify the size of the waterbox explicitly or calculate it 611 automatically for a macromolecule along with adding padding 612 around ther macromolecule. 613 size: A space delimited triplet of values corresponding to water 614 size in nanometers. It must be specified during 'Size' value of 615 'mode' parameter. 616 padding: Padding around a macromolecule in nanometers for filling 617 box with water. It must be specified during 'Padding' value of 618 'mode' parameter. 619 620 -w, --workingdir <dir> 621 Location of working directory which defaults to the current directory. 622 623 Examples: 624 To prepare a macromolecule in a PDB file by replacing non-standard residues, 625 adding missing residues, adding missing heavy atoms and missing hydrogens, 626 and generate a PDB file, type: 627 628 % OpenMMPrepareMacromolecule.py -i Sample11.pdb -o Sample11Out.pdb 629 630 To run the first example for listing additional details about missing atoms and 631 residues, and generate a PDB file, type: 632 633 % OpenMMPrepareMacromolecule.py --listDetails yes -i Sample11.pdb 634 -o Sample11Out.pdb 635 636 To run the first example for deleting all heterogens including water along 637 with performing all default actions, and generate a PDB file, type: 638 639 % OpenMMPrepareMacromolecule.py --deleteHeterogens All -i Sample11.pdb 640 -o Sample11Out.pdb 641 642 To run the first example for deleting water only heterogens along with performing 643 all default actions, and generate a PDB file, type: 644 645 % OpenMMPrepareMacromolecule.py --deleteHeterogens WaterOnly 646 -i Sample11.pdb -o Sample11Out.pdb --ov 647 648 To run the first example for adding a water box by automatically calculating its 649 size, along with performing all default actions, and generate a PDB file, type: 650 651 % OpenMMPrepareMacromolecule.py --waterBox yes -i Sample11.pdb 652 -o Sample11Out.pdb 653 654 To run the previous example by explcitly specifying various water box parameters, 655 and generate a PDB file, type: 656 657 % OpenMMPrepareMacromolecule.py --waterBox yes 658 --waterBoxParams "mode,Padding, padding, 1.0, shape, cube" 659 -i Sample11.pdb -o Sample11Out.pdb 660 661 To run the first example for adding a water box of a specific size along with 662 performing all default actions, and generate a PDB file, type: 663 664 % OpenMMPrepareMacromolecule.py --waterBox yes 665 --waterBoxParams "mode,Size, size, 7.635 7.077 7.447" 666 -i Sample11.pdb -o Sample11Out.pdb 667 668 To add a lipid membrane around a correctly oriented and positioned macromolecule, 669 along with performing all default actions, and generate a PDB file, type: 670 671 % OpenMMPrepareMacromolecule.py --membrane yes 672 -i Sample12.pdb -o Sample12Out.pdb 673 674 To add a lipid membrane around a correctly oriented and positioned macromolecule, 675 deleting all heterogens along with performing all default actions, and generate a 676 PDB file, type: 677 678 % OpenMMPrepareMacromolecule.py --membrane yes --deleteHeterogens All 679 -i Sample12.pdb -o Sample12Out.pdb 680 681 To run the previous example by explcitly specifying various membrane parameters, 682 and generate a PDB file, type: 683 684 % OpenMMPrepareMacromolecule.py --membrane yes 685 --membraneParams "lipidType, POPC, membraneCenterZ, 0.0, padding, 1.0" 686 -i Sample12.pdb -o Sample12Out.pdb 687 688 Author: 689 Manish Sud(msud@san.rr.com) 690 691 See also: 692 InfoPDBFiles.pl, ExtractFromPDBFiles.pl, PyMOLExtractSelection.py, 693 PyMOLInfoMacromolecules.py, PyMOLSplitChainsAndLigands.py 694 695 Copyright: 696 Copyright (C) 2025 Manish Sud. All rights reserved. 697 698 The functionality available in this script is implemented using OpenMM, an 699 open source molecuar simulation package. 700 701 This file is part of MayaChemTools. 702 703 MayaChemTools is free software; you can redistribute it and/or modify it under 704 the terms of the GNU Lesser General Public License as published by the Free 705 Software Foundation; either version 3 of the License, or (at your option) any 706 later version. 707 708 """ 709 710 if __name__ == "__main__": 711 main()
