Reorientation of glycoprotein for umbrella sampling and where the imaginary spring is attached

[Johnla]

GROMACS version:2021.6
GROMACS modification: Yes/No
Here post your question
Hello guys, I am trying to do an umbrella sampling by pulling the binding partner from the target protein using their complex structure.
My first question is how can I create an ideal box for umbrella sampling? As the target protein has extensive glyosylations I need to use CHARMM-GUI to prepare the input for gromacs. When it was finished, I found the orientation of the complex was not ideal for pulling as the imaginary pulling coordinate was not parallele to one of there axes. I could imagine that the binding partner would quickly move out of the box during the pull. So do I have to manually ajdust the orientation of the complex using either pymol or vmd before submitting to the CHARMM-GUI or there’s some other tools that can help doing it automatically and more precisely?

My second question is about the pulling itself. Through learning Justin’s umbrella sampling tutorial. I kind of get the idea of specifying the two groups “Chain_A” and “Chain_B”, and I konw “Chain_B” is restrained by a static umbrella forces. So I want to know how the progrom determine where the imaginary spring is attached? It looks to me that once you specify “Chain_B” to be restrained, “Chain_A” was automatically chosen by the program to attach the imaginary spring. So do I understand this right?

My third question is about where on “Chain_A” the imaginary spring is attached. In the tutorial, I didn’t find the parameter specifying this. The cartoon in Justin’s tutorial looks as if the imaginary spring is attached to the terminal of “Chain_A”. When I used vmd to visualize this, it gave me a similar video like the one Justin put on the web. But when I used pymol to visualize this, it gave me a different video like the "imaginary spring was attached to the center of the “Chain_A”. How could this happen?

[scinikhil]

You need to orient the molecule or system based on the appropriate direction you wanted to pull. you are on the right direction by ortienting pymol or any packages that can do that. when you do that you make sure that the ‘newstrucrure’ coordinates are the one is saved.

The concept of imaginary spring is in the mdp file, Chain_A is defined so, if not specified the center of mass of the specific group is considered. If you wanted to get a better view, you may specify the pull-vector or group in the mdp file. You need to carefully define the vector or the direction to see the desired movement of the group, also the spring constant is another way to control the movement.

[Johnla]

Hi scinikhil,
I tried CHARMM-GUI. But still there’s another problem. In ChARMM-GUI there’s no place to specify where the complex should be placed in the box. I just run through the “Glycan Reader and Modeller” and it finally generated a .gro file that placed the complex in the center of the box. In gromacs we can use “gmx editconf” to specify the center of the complex. But in the umbrella sampling tutorial, this is done before puting water and ion molecules. So do you have any suggestions to solve this problem?
Thank you!

[Johnla]

Hello guys,
I’ve figured out the correct way to do this.
1, Orient the pdb file in pymol using the “rotate” command to rotate the complex in the suitable orientation for pulling. E.g., “roate x, sel, 90” will rotate the selected chain around the x direction 90 degree. After inspecting mine again carefully I found the x axis is suitable for pulling, so I didn’t rotate the complex.
2, In CHARMM-GUI, select “Input generator” and then “Glycan reader/modeler”, Upload the structure to the CHARMM-GUI. The program willl automactically detect the glycan chains in the protein. Check if the glycan were correctly represented. If not, there might be a TER between the protein chain and the glycan, delete this TER will correct this problem.
3. Just follow the program sequence and press NEXT STEP. Under the “PDB info” tab, make sure the correct disulfide bond linkage are selected if the proteins contains disufide bonds.
4. Under the “CHARMM-PDB” tab, check “Fit water box size”, and “Enter the edge distance” automactically is filled with 10, which is by default, the unit I believe is angstrom. In this way, a box size is automatically chosen according to the size of the complex which will be placed at the center of box. Filled in desired ions and concentrations.
5. Make sure “gromacs” as the “the output format was chosen.
6. Download the generated files as a compressed file. Uncompress the file using the “tar xvzf filename.gz”
7. cd into the uncompressed folder, you will find a folder named gromacs, cd into that folder and there’s a file named “step3_input.gro”, copy this file and topol.top files into your working folder.
8. Using a texeditor to open the step3_input.gro file, delete any line containing the TIP3, CLA, and SOD or POT, count the rest of the lines number using “cat filename.gro | wc -l”. Coordinate_number = line_number - 3, Edit the .gro file again change the total atom number to the Coordinate_number. Save it as another file say, “input_file.gro” and take a note of the box size info
9. Now, using the gmx editconf command to make a new box of the desired dimention. I choose the same y and z length as the one generated by CHARMM-GUI, but doubled the x length, which I think might be enough to pull the molecule.

10. Using the gmx solvate command to fill water into this box.
11. Using the gmx make_ndx command to set the index name for the two component of the complex, say “PROA” and PROB”. gmx make_ndx -f solved.gro Attention should be taken here that the residue Number will be set individually after “gmx genion”, making selection of residue in making the index confusing. To correctly select the residue number, it is recommended to use the solved.gro file from the output file of the gmx solvate command or the output .gro. file from the CHARMM-GUI.
12. Using the command line to calculate the number of solvant molecules “grep SOL solved.gro | wc -l”. Divide the output number by 3 give the solvant molecule number. Using a text editor to open the topol.top file, delete the ions such as “CLA”, “SOD”, “POT” and change the TIP3 number to the calculaed solvant number in the SOL solved.gro file.
under the [molecule], save this file as “solved_corrected.gro”.
from:

to:

Save the topol.top file.
13. Copy the “toppar” folder in the CHARMM-GUI generated in the working folder. Use the gmx grompp to prepare the tpr file for gmx genion to replace the solvant molecule (“SOL”) in the solved.gro file with “TIP3” Because the solvant name in the .gro file is different as in the index file, which is TIP3, SOL

When running the gmx grompp for ion preparation, the system will report warnings and fail this step, add the “maxwarn -1” option to negect the warnings. I am using gromacs 2021.6, newer versions of gromacs may use a large number for this option, such as “maxwarn 400”
gmx grompp -f ions.mdp -c aCTX_a1ECD_noglycan_modeling_solved.gro -p topol.top -o ions.tpr -maxwarn -1
14. Run "gmx genion, making sure the same ion name is used as in the tpr files generated by CHARMM-GUI in the toppar folder. In Justin’s tutorial, NA, CL were used for sodium and chloride, but the CHARMM-GUI use SOD and CLA, respectively.
gmx genion -s ions.tpr -o solved_ion.gro -p topol.top -pname SOD -nname CLA -neutral -conc 0.15
15. Run energy minimization using the em.mdp file in Justin’s tutorial
16. Edit the Energy restraint in the proteins topology file in the toppar folder, replace the POSRES_FC_BB, POSRES_FC_SC to 1000
from:

to:

17. Run equilibration using the nvt.mdp and npt.mdp file provided in Justin’s lysozyme tutorial.

gmx make_ndx -f solved.gro
18. Edit the md_pull.mdp file, change the index name of the two groups involved in the pulling simulation, also make sure the correct pulling axis is selected

I will update this thread when I make further progress.