Adapt protein mutation setup pipeline to handle glycosylated RBD:ACE2

perses/data/rbd-ace2/1_rbd_ace2_template_tleap.in

@@ -0,0 +1,131 @@
+source leaprc.GLYCAM_06j-1 #  Source leaprc file for GLYCAM force field
+source leaprc.protein.ff14SB #  Source leaprc file for ff14SB protein force field
+source leaprc.water.tip3p #  Source leaprc file for TIP3P water model (and ions)
+
+
+######### RBD #########
+
+mol1 = loadpdb 0_rbd_tleap.pdb
+
+#  Add disulphide bridges
+
+bond mol1.336.SG mol1.361.SG
+bond mol1.379.SG mol1.432.SG
+bond mol1.391.SG mol1.525.SG
+bond mol1.480.SG mol1.488.SG
+
+#  N343 glycans
+
+bond mol1.343.ND2 mol1.528.C1 #  Bond N343 to GlcNAc
+bond mol1.528.O6 mol1.537.C1
+bond mol1.528.O4 mol1.529.C1
+bond mol1.529.O4 mol1.530.C1
+bond mol1.530.O6 mol1.531.C1
+bond mol1.531.O2 mol1.532.C1
+bond mol1.532.O4 mol1.533.C1
+bond mol1.530.O3 mol1.534.C1
+bond mol1.534.O2 mol1.535.C1
+bond mol1.535.O4 mol1.536.C1
+
+######### ACE2 #########
+
+mol2 = loadpdb 0_ace2_tleap.pdb
+
+#  Add disulphide bridges
+bond mol2.133.SG mol2.141.SG
+bond mol2.344.SG mol2.361.SG
+bond mol2.530.SG mol2.542.SG
+
+#  N53 glycans
+
+bond mol2.53.ND2 mol2.727.C1 #  Bond ASN53 to GlcNAc
+bond mol2.727.O6 mol2.734.C1
+bond mol2.727.O4 mol2.728.C1
+bond mol2.728.O4 mol2.729.C1
+bond mol2.729.O6 mol2.730.C1
+bond mol2.730.O2 mol2.731.C1
+bond mol2.729.O3 mol2.732.C1
+bond mol2.732.O2 mol2.733.C1
+
+#  N90 glycans
+
+bond mol2.90.ND2 mol2.735.C1 #  Bond ASN90 to GlcNAc
+bond mol2.735.O6 mol2.743.C1
+bond mol2.735.O4 mol2.736.C1
+bond mol2.736.O4 mol2.737.C1
+bond mol2.737.O6 mol2.738.C1
+bond mol2.738.O2 mol2.739.C1
+bond mol2.739.O4 mol2.740.C1
+bond mol2.737.O3 mol2.741.C1
+bond mol2.741.O2 mol2.742.C1
+
+#  N103 glycans
+
+bond mol2.103.ND2 mol2.744.C1 #  Bond ASN103 to GlcNAc
+bond mol2.744.O6 mol2.751.C1
+bond mol2.744.O4 mol2.745.C1
+bond mol2.745.O4 mol2.746.C1
+bond mol2.746.O6 mol2.747.C1
+bond mol2.747.O2 mol2.748.C1
+bond mol2.746.O3 mol2.749.C1
+bond mol2.749.O2 mol2.750.C1
+
+#  N322 glycans
+
+bond mol2.322.ND2 mol2.752.C1 #  Bond ASN322 to GlcNAc
+bond mol2.752.O6 mol2.759.C1
+bond mol2.752.O4 mol2.753.C1
+bond mol2.753.O4 mol2.754.C1
+bond mol2.754.O6 mol2.755.C1
+bond mol2.755.O2 mol2.756.C1
+bond mol2.754.O3 mol2.757.C1
+bond mol2.757.O2 mol2.758.C1
+
+# N432 glycans
+
+bond mol2.432.ND2 mol2.760.C1 #  Bond ASN432 to GlcNAc
+bond mol2.760.O6 mol2.769.C1
+bond mol2.760.O4 mol2.761.C1
+bond mol2.761.O4 mol2.762.C1
+bond mol2.762.O6 mol2.763.C1
+bond mol2.763.O2 mol2.764.C1
+bond mol2.764.O4 mol2.765.C1
+bond mol2.762.O3 mol2.766.C1
+bond mol2.766.O2 mol2.767.C1
+bond mol2.767.O4 mol2.768.C1
+
+# N546 glycans
+
+bond mol2.546.ND2 mol2.770.C1 #  Bond ASN546 to GlcNAc
+bond mol2.770.O4 mol2.771.C1
+bond mol2.771.O4 mol2.772.C1
+bond mol2.772.O6 mol2.773.C1
+bond mol2.773.O2 mol2.774.C1
+bond mol2.772.O3 mol2.775.C1
+bond mol2.775.O2 mol2.776.C1
+
+# N690 glycans
+
+bond mol2.690.ND2 mol2.777.C1 #  Bond ASN690 to GlcNAc
+bond mol2.777.O6 mol2.784.C1
+bond mol2.777.O4 mol2.778.C1
+bond mol2.778.O4 mol2.779.C1
+bond mol2.779.O6 mol2.780.C1
+bond mol2.780.O2 mol2.781.C1
+bond mol2.779.O3 mol2.782.C1
+bond mol2.782.O2 mol2.783.C1
+
+# Combine both systems to make RBD:ACE2 complex
+complex = combine {mol1 mol2}
+
+# Solvate the system and add ions
+charge complex
+solvateoct complex TIP3PBOX 9 #  Add some solvent with padding = 9
+addionsrand complex Na+ 187 Cl- 167 #  Randomly replace waters with ions, waters = 61736, see README for how No. ions chosen
+
+# Save everything
+savepdb complex 1_rbd_ace2_tleap.pdb
+
+saveamberparm complex 1_rbd_ace2_tleap.prmtop 2_rbd_ace2_tleap.inpcrd
+
+quit #  Quit tleap program

perses/data/rbd-ace2/1_rbd_template_tleap.in

@@ -0,0 +1,40 @@
+source leaprc.GLYCAM_06j-1 #  Source leaprc file for GLYCAM force field
+source leaprc.protein.ff14SB #  Source leaprc file for ff14SB protein force field
+source leaprc.water.tip3p #  Source leaprc file for TIP3P water model (and ions)
+
+
+######### RBD #########
+
+mol1 = loadpdb 0_rbd_tleap.pdb
+
+#  Add disulphide bridges
+
+bond mol1.336.SG mol1.361.SG
+bond mol1.379.SG mol1.432.SG
+bond mol1.391.SG mol1.525.SG
+bond mol1.480.SG mol1.488.SG
+
+#  N343 glycans
+
+bond mol1.343.ND2 mol1.528.C1 #  Bond N343 to GlcNAc
+bond mol1.528.O6 mol1.537.C1
+bond mol1.528.O4 mol1.529.C1
+bond mol1.529.O4 mol1.530.C1
+bond mol1.530.O6 mol1.531.C1
+bond mol1.531.O2 mol1.532.C1
+bond mol1.532.O4 mol1.533.C1
+bond mol1.530.O3 mol1.534.C1
+bond mol1.534.O2 mol1.535.C1
+bond mol1.535.O4 mol1.536.C1
+
+# Solvate the system and add ions
+charge mol1
+solvateoct mol1 TIP3PBOX 9 #  Add some solvent with padding = 9
+addionsrand mol1 Na+ 50 Cl- 52 #  Randomly replace waters with ions, waters = 18534, see README for how No. ions chosen
+
+# Save everything
+savepdb mol1 0_rbd_tleap.pdb
+
+saveamberparm mol1 0_rbd_tleap.prmtop 0_rbd_tleap.inpcrd
+
+quit #  Quit tleap program

perses/data/rbd-ace2/3_mutate.py

@@ -0,0 +1,51 @@
+import pymol
+from pymol import cmd
+import sys
+import argparse
+
+# https://sourceforge.net/p/pymol/mailman/message/27979284/
+# https://sourceforge.net/p/pymol/mailman/message/11671708/
+
+# Read in arguments
+parser = argparse.ArgumentParser(description='run pymol mutagenesis')
+parser.add_argument('input_pdb', type=str, help='path to input file')
+parser.add_argument('output_pdb', type=str, help='path to output file')
+parser.add_argument('selection', type=str, help='e.g. "R/337/"')
+parser.add_argument('mutant', type=str, help='amino acid three letter code, e.g. ALA')
+args = parser.parse_args()
+
+# Beware, the error message you get is because pymol unsuccessfully tries 
+#to open the command line arguments as files.
+
+d = {'CYS': 'C', 'ASP': 'D', 'SER': 'S', 'GLN': 'Q', 'LYS': 'K',
+     'ILE': 'I', 'PRO': 'P', 'THR': 'T', 'PHE': 'F', 'ASN': 'N', 
+     'GLY': 'G', 'HIS': 'H', 'LEU': 'L', 'ARG': 'R', 'TRP': 'W', 
+     'ALA': 'A', 'VAL':'V', 'GLU': 'E', 'TYR': 'Y', 'MET': 'M'}
+
+print(f"Mutating {args.mutant}...")
+
+# Launch pymol session
+pymol.pymol_argv = ["pymol", "-qc"] + sys.argv[1:]
+pymol.finish_launching()
+
+# Load RBD (no solvent)
+cmd.load(args.input_pdb)
+
+# Mutate
+cmd.wizard("mutagenesis")
+cmd.do("refresh_wizard")
+cmd.get_wizard().set_mode(args.mutant)
+cmd.get_wizard().do_select(args.selection) 
+
+# Select rotamer
+cmd.frame(1)
+
+# Apply the mutation
+cmd.get_wizard().apply()
+cmd.set_wizard() # Equivalent to clicking "Done" in the GUI
+
+# Save
+cmd.save(args.output_pdb)
+# Wait for everything to be done
+cmd.refresh() 
+

perses/data/rbd-ace2/5_rbd_ace2_mutant_template_tleap.in

@@ -0,0 +1,128 @@
+source leaprc.GLYCAM_06j-1 #  Source leaprc file for GLYCAM force field
+source leaprc.protein.ff14SB #  Source leaprc file for ff14SB protein force field
+source leaprc.water.tip3p #  Source leaprc file for TIP3P water model (and ions)
+
+######### RBD #########
+
+mol1 = loadpdb 4_rbd_ace2_mutant_rbd_tleap_final.pdb
+#  Add disulphide bridges
+
+bond mol1.336.SG mol1.361.SG
+bond mol1.379.SG mol1.432.SG
+bond mol1.391.SG mol1.525.SG
+bond mol1.480.SG mol1.488.SG
+
+#  N343 glycans
+
+bond mol1.343.ND2 mol1.528.C1 #  Bond N343 to GlcNAc
+bond mol1.528.O6 mol1.537.C1
+bond mol1.528.O4 mol1.529.C1
+bond mol1.529.O4 mol1.530.C1
+bond mol1.530.O6 mol1.531.C1
+bond mol1.531.O2 mol1.532.C1
+bond mol1.532.O4 mol1.533.C1
+bond mol1.530.O3 mol1.534.C1
+bond mol1.534.O2 mol1.535.C1
+bond mol1.535.O4 mol1.536.C1
+
+######### ACE2 #########
+
+mol2 = loadpdb 4_rbd_ace2_mutant_ace2_tleap_final.pdb
+
+#  Add disulphide bridges
+bond mol2.133.SG mol2.141.SG
+bond mol2.344.SG mol2.361.SG
+bond mol2.530.SG mol2.542.SG
+
+#  N53 glycans
+
+bond mol2.53.ND2 mol2.727.C1 #  Bond ASN53 to GlcNAc
+bond mol2.727.O6 mol2.734.C1
+bond mol2.727.O4 mol2.728.C1
+bond mol2.728.O4 mol2.729.C1
+bond mol2.729.O6 mol2.730.C1
+bond mol2.730.O2 mol2.731.C1
+bond mol2.729.O3 mol2.732.C1
+bond mol2.732.O2 mol2.733.C1
+
+#  N90 glycans
+
+bond mol2.90.ND2 mol2.735.C1 #  Bond ASN90 to GlcNAc
+bond mol2.735.O6 mol2.743.C1
+bond mol2.735.O4 mol2.736.C1
+bond mol2.736.O4 mol2.737.C1
+bond mol2.737.O6 mol2.738.C1
+bond mol2.738.O2 mol2.739.C1
+bond mol2.739.O4 mol2.740.C1
+bond mol2.737.O3 mol2.741.C1
+bond mol2.741.O2 mol2.742.C1
+
+#  N103 glycans
+
+bond mol2.103.ND2 mol2.744.C1 #  Bond ASN103 to GlcNAc
+bond mol2.744.O6 mol2.751.C1
+bond mol2.744.O4 mol2.745.C1
+bond mol2.745.O4 mol2.746.C1
+bond mol2.746.O6 mol2.747.C1
+bond mol2.747.O2 mol2.748.C1
+bond mol2.746.O3 mol2.749.C1
+bond mol2.749.O2 mol2.750.C1
+
+#  N322 glycans
+
+bond mol2.322.ND2 mol2.752.C1 #  Bond ASN322 to GlcNAc
+bond mol2.752.O6 mol2.759.C1
+bond mol2.752.O4 mol2.753.C1
+bond mol2.753.O4 mol2.754.C1
+bond mol2.754.O6 mol2.755.C1
+bond mol2.755.O2 mol2.756.C1
+bond mol2.754.O3 mol2.757.C1
+bond mol2.757.O2 mol2.758.C1
+
+# N432 glycans
+
+bond mol2.432.ND2 mol2.760.C1 #  Bond ASN432 to GlcNAc
+bond mol2.760.O6 mol2.769.C1
+bond mol2.760.O4 mol2.761.C1
+bond mol2.761.O4 mol2.762.C1
+bond mol2.762.O6 mol2.763.C1
+bond mol2.763.O2 mol2.764.C1
+bond mol2.764.O4 mol2.765.C1
+bond mol2.762.O3 mol2.766.C1
+bond mol2.766.O2 mol2.767.C1
+bond mol2.767.O4 mol2.768.C1
+
+# N546 glycans
+
+bond mol2.546.ND2 mol2.770.C1 #  Bond ASN546 to GlcNAc
+bond mol2.770.O4 mol2.771.C1
+bond mol2.771.O4 mol2.772.C1
+bond mol2.772.O6 mol2.773.C1
+bond mol2.773.O2 mol2.774.C1
+bond mol2.772.O3 mol2.775.C1
+bond mol2.775.O2 mol2.776.C1
+
+# N690 glycans
+
+bond mol2.690.ND2 mol2.777.C1 #  Bond ASN690 to GlcNAc
+bond mol2.777.O6 mol2.784.C1
+bond mol2.777.O4 mol2.778.C1
+bond mol2.778.O4 mol2.779.C1
+bond mol2.779.O6 mol2.780.C1
+bond mol2.780.O2 mol2.781.C1
+bond mol2.779.O3 mol2.782.C1
+bond mol2.782.O2 mol2.783.C1
+
+######### Solvent #########
+mol3 = loadpdb 4_rbd_ace2_mutant_solvent_tleap_final.pdb
+
+# Combine both molecules to make solvated RBD:ACE2 complex
+complex = combine {mol1 mol2 mol3}
+
+# Add box dimensions
+set complex box {114.609 118.869 170.037}
+
+# Save everything
+savepdb complex 5_rbd_ace2_mutant_tleap.pdb
+saveamberparm complex 5_rbd_ace2_mutant_tleap.prmtop 5_rbd_ace2_mutant_tleap.inpcrd
+quit #  Quit tleap program

perses/data/rbd-ace2/5_rbd_mutant_template_tleap.in

@@ -0,0 +1,41 @@
+source leaprc.GLYCAM_06j-1 #  Source leaprc file for GLYCAM force field
+source leaprc.protein.ff14SB #  Source leaprc file for ff14SB protein force field
+source leaprc.water.tip3p #  Source leaprc file for TIP3P water model (and ions)
+
+
+######### RBD #########
+
+mol1 = loadpdb 4_rbd_mutant_rbd_tleap_final.pdb
+#  Add disulphide bridges
+
+bond mol1.336.SG mol1.361.SG
+bond mol1.379.SG mol1.432.SG
+bond mol1.391.SG mol1.525.SG
+bond mol1.480.SG mol1.488.SG
+
+#  N343 glycans
+
+bond mol1.343.ND2 mol1.528.C1 #  Bond N343 to GlcNAc
+bond mol1.528.O6 mol1.537.C1
+bond mol1.528.O4 mol1.529.C1
+bond mol1.529.O4 mol1.530.C1
+bond mol1.530.O6 mol1.531.C1
+bond mol1.531.O2 mol1.532.C1
+bond mol1.532.O4 mol1.533.C1
+bond mol1.530.O3 mol1.534.C1
+bond mol1.534.O2 mol1.535.C1
+bond mol1.535.O4 mol1.536.C1
+
+######### Solvent #########
+mol2 = loadpdb 4_rbd_mutant_solvent_tleap_final.pdb
+
+# Combine both molecules to make solvated RBD
+complex = combine {mol1 mol2}
+
+# Add box dimensions
+set complex box {114.609 118.869 170.037}
+
+# Save everything
+savepdb complex 5_rbd_mutant_tleap.pdb
+saveamberparm complex 5_rbd_mutant_tleap.prmtop 5_rbd_mutant_tleap.inpcrd
+quit #  Quit tleap program