use yaml file for metadynamics params

enhanced_sampling/cv_building.py

@@ -1,5 +1,4 @@
-import openmm, openmmtools
-import mdtraj
+import openmm
 
 # def create_rmsd_restraint(positions, atom_indicies):
 #     rmsd_cv = openmm.RMSDForce(positions, atom_indicies)

hpc_scripts/rmsd_metadynamics.sh

@@ -20,4 +20,5 @@
 
 source ~/.bashrc
 conda activate enhanced-sampling-tmem175
-python ../metadynamics/rmsd_force.py -i ../systems/system00 -r ../systems/system01 -o ../../enhanced-sampling-tmem175-data/rmsd_metadynamics
+python ../metadynamics/rmsd_force.py -i ../systems/system00 -r ../systems/system01 -o ../../enhanced-sampling-tmem175-data/rmsd_metadynamics -y ../metadynamics/tic1_0.7.yaml
+#python ../metadynamics/rmsd_force.py -i ../systems/system00 -r ../systems/system01 -o ../../enhanced-sampling-tmem175-data/rmsd_metadynamics -y ../metadynamics/cas.yaml

log/2022-09-13_log.yaml

@@ -108,3 +108,51 @@
   notes: ''
   script: rmsd_force.py
   success: ''
+'2022-09-13T16:59:24.717987':
+  args:
+    charmm_param_dir: ../forcefield_data/charmm-params
+    input_dir: ../systems/system00
+    output_dir: /Volumes/Rohirrim/local_test/enhanced_sampling/metadynamics_rmsd
+    params_file: ../sim_params/local_test.yaml
+    reference_dir: ../systems/system01
+    res_list: /Users/alexpayne/Scientific_Projects/fah-analysis-tmem175/metadata/key_res.yaml
+    rmsd_sel: protein_ca
+  notes: ''
+  script: rmsd_force.py
+  success: ''
+'2022-09-13T18:27:19.218283':
+  args:
+    charmm_param_dir: ../forcefield_data/charmm-params
+    input_dir: ../systems/system00
+    meta_params: tic1_0.7.yaml
+    output_dir: /Volumes/Rohirrim/local_test/enhanced_sampling/metadynamics_rmsd
+    params_file: ../sim_params/local_test.yaml
+    reference_dir: ../systems/system01
+    rmsd_sel: protein_ca
+  notes: ''
+  script: rmsd_force.py
+  success: ''
+'2022-09-13T18:31:35.737518':
+  args:
+    charmm_param_dir: ../forcefield_data/charmm-params
+    input_dir: ../systems/system00
+    meta_params: tic1_0.7.yaml
+    output_dir: /Volumes/Rohirrim/local_test/enhanced_sampling/metadynamics_rmsd
+    params_file: ../sim_params/local_test.yaml
+    reference_dir: ../systems/system01
+    rmsd_sel: protein_ca
+  notes: ''
+  script: rmsd_force.py
+  success: ''
+'2022-09-13T18:33:28.228150':
+  args:
+    charmm_param_dir: ../forcefield_data/charmm-params
+    input_dir: ../systems/system00
+    meta_params: cas.yaml
+    output_dir: /Volumes/Rohirrim/local_test/enhanced_sampling/metadynamics_rmsd
+    params_file: ../sim_params/local_test.yaml
+    reference_dir: ../systems/system01
+    rmsd_sel: protein_ca
+  notes: ''
+  script: rmsd_force.py
+  success: ''

metadynamics/cas.yaml

@@ -0,0 +1,4 @@
+selection: protein_ca
+min_value: 0.1
+max_value: 3.0
+bias_width: 0.02

metadynamics/rmsd_force.py

@@ -5,9 +5,11 @@
 import openmm.app
 import logging
 
+import yaml
+
 repo_path = f"{os.path.dirname(os.path.dirname(os.path.abspath(__file__)))}"
 sys.path.append(repo_path)
-from enhanced_sampling import utils, system_building as sb, system_saving as ss, reporters
+from enhanced_sampling import utils, system_building as sb, system_saving as ss, reporters, cv_building
 
 
 ## ARGUMENT PARSING
@@ -23,6 +25,8 @@ def get_args():
                         help="Directory containing charmm parameter files")
     parser.add_argument("-p", "--params_file", type=str, default="../sim_params/defaults.yaml",
                         help="YAML file containing simulation running values")
+    parser.add_argument('-y', "--meta_params", type=str, default=None,
+                        help="Path to yaml file containing residue list and min_max values")
     args = parser.parse_args()
     return args
 
@@ -83,15 +87,28 @@ def main():
 
     assert len(ref_positions) == len(input_dict["positions"])
 
-    idx = [atom.index for atom in ref_psf.topology.atoms() if
-           atom.residue.chain.index == 0 or atom.residue.chain.index == 1]
+    if args.meta_params:
+        assert os.path.exists(args.meta_params)
+        with open(args.meta_params, "r") as f:
+            meta_dict = yaml.safe_load(f)
+        res_list = meta_dict.get('res_list')
+        rmsd_sel = meta_dict.get('selection')
+        min_value = meta_dict.get('min_value')
+        max_value = meta_dict.get('max_value')
+        bias_width = meta_dict.get('bias_width')
+    else:
+        min_value = 0.01
+        max_value = 0.4
+        rmsd_sel = 'protein_ca'
+
+    idx = cv_building.get_openmm_idx(ref_psf.topology, rmsd_sel, res_list)
 
     rmsd_force = openmm.RMSDForce(ref_positions, idx)
 
     rmsd_bias = openmm.app.metadynamics.BiasVariable(force=rmsd_force,
-                                                     minValue=0.01,
-                                                     maxValue=0.4,
-                                                     biasWidth=0.02,
+                                                     minValue=min_value,
+                                                     maxValue=max_value,
+                                                     biasWidth=bias_width,
                                                      periodic=False)
 
     meta = openmm.app.Metadynamics(system=system,

metadynamics/tic1_0.7.yaml

@@ -0,0 +1,58 @@
+selection: protein_heavy
+min_value: 0.2
+max_value: 8.0
+bias_width: 0.05
+res_list:
+  - 46
+  - 55
+  - 56
+  - 57
+  - 58
+  - 59
+  - 60
+  - 61
+  - 65
+  - 66
+  - 70
+  - 71
+  - 74
+  - 267
+  - 271
+  - 273
+  - 274
+  - 275
+  - 276
+  - 278
+  - 279
+  - 280
+  - 283
+  - 308
+  - 309
+  - 312
+  - 316
+  - 317
+  - 320
+  - 356
+  - 357
+  - 358
+  - 359
+  - 360
+  - 361
+  - 362
+  - 363
+  - 364
+  - 365
+  - 366
+  - 367
+  - 368
+  - 369
+  - 370
+  - 371
+  - 372
+  - 373
+  - 375
+  - 376
+  - 440
+  - 441
+  - 442
+  - 443