Frequently Asked Questions¶
I installed the software, how do I test if it is correctly installed?¶
Because different users have different settings and requirements for their clusters or workstations, we provide a general job handling script for you to customize to your needs.
To assist with job handling script customization, example input files are
available under the $SILCSBIODIR/examples folder.
For SILCS, use the following commands to make sure the software is correctly installed and the job handling script is working. If you are only interested in SSFEP simulations, you may skip to the SSFEP section below.
mkdir -p test/silcs
cd test/silcs
cp $SILCSBIODIR/examples/silcs/p38a.pdb .
$SILCSBIODIR/silcs/1_setup_silcs_boxes prot=p38a.pdb
$SILCSBIODIR/silcs/2a_run_gcmd prot=p38a.pdb numsys=1 nproc=1
If this set of commands runs without error, confirm that the SILCS job is
running with the check_progress command:
$SILCSBIODIR/silcs/check_progress
and then go ahead and stop the successfully running SILCS job:
$SILCSBIODIR/silcs/2a_run_gcmd cancel=true sys=1
and confirm it is stopped:
$SILCSBIODIR/silcs/check_progress
Otherwise, if you experienced an error with the 1_setup-silcs_boxes step,
the software is not correctly installed, whereas if you experienced an error
with the 2a_run_gcmd step, the job handling scripts need to be edited. The
job handling scripts for SILCS are:
templates/silcs/job_mc_md.tmpltemplates/silcs/job_gen_maps.tmpltemplates/silcs/pymol_fragmap.tmpltemplates/silcs/vmd_fragmap.tmpltemplates/silcs/job_cleanup.tmpl
Typically the header portion of a job handling script requires editing. Please contact support@silcsbio.com if you need assistance.
For SSFEP, use the following commands to make sure the software is correctly installed and the job handling script is working.
mkdir -p test/ssfep
cd test/ssfep
cp $SILCSBIODIR/examples/ssfep/* .
$SILCSBIODIR/ssfep/1_setup_ssfep prot=4ykr.pdb lig=lig.mol2
$SILCSBIODIR/silcs/2_run_md_ssfep prot=4ykr.pdb lig=lig.mol2 nproc=1
If this set of commands runs without error, confirm that the SSFEP job is
running with the check_progress command:
$SILCSBIODIR/ssfep/check_progress
and then go ahead and stop the successfully running SSFEP job:
$SILCSBIODIR/ssfep/2_run_md_ssfep cancel=true target=lig
$SILCSBIODIR/ssfep/2_run_md_ssfep cancel=true target=prot
and confirm it is stopped:
$SILCSBIODIR/ssfep/check_progress
Otherwise, if you experienced an error with the 1_setup_ssfep step, the
software is not correctly installed, and if you experienced an error with the
2_run_md_ssfep step, the job handling script needs to be edited. The job
handling scripts for SSFEP are:
templates/ssfep/job_lig_md.tmpltemplates/ssfep/job_prot_lig_md.tmpltemplates/ssfep/job_dG.tmpl
Typically the header portion of a job handling script requires editing. Please contact support@silcsbio.com if you need assistance.
I don’t have a cluster but I have a GPU workstation. What can I do?¶
You may be able to practically run the SilcsBio software if your GPU workstation has sufficient resources. An appropriate workstation may have at least 24 CPU cores, 4 GPUs, 64 GB of RAM, and 10 TB of disk space. Installing a job queueing system, such as the Slurm Workload Manager, will allow the SilcsBio server software to run on the workstation.
The SilcsBio Workstation is a turn-key GPU workstation hardware+software solution developed by SilcsBio that comes with all necessary software pre-installed. The SilcsBio workstation has a quiet, sleek form factor for use in an office setting and comes ready to plug in to a standard wall electrical socket. Please contact info@silcsbio.com for details.
I compiled my GROMACS with MPI and my job is not running¶
Please contact us so we can repackage the files with the appropriate
command using mpirun instead.
Alternatively, you may edit the job handling script to edit the GROMACS command.
For example, the mdrun command is specified at the top of
templates/ssfep/job_lig_md.tmpl file:
mdrun="${GMXDIR}/gmx mdrun -nt $nproc"
You may edit this to
mdrun="mpirun -np $nproc ${GMXDIR}/gmx mdrun"
GROMACS on the head node does not run because the head node and compute node have different operating systems¶
In this case, we recommend compilng GROMACS on the head node and
compiling mdrun only on the compute node.
Building only mdrun can be done by supplying the
-DGMX_BUILD_MDRUN_ONLY=on keyword to the cmake command in the build
process. Once the mdrun program is built, place it in the same $GMXDIR
folder. Now template files needs to be edited to use the mdrun command
properly on the compute node.
For example, the mdrun command is specified at the top of the
templates/ssfep/job_lig_md.tmpl file:
mdrun="${GMXDIR}/gmx mdrun -nt $nproc"
You may edit this to
mdrun="${GMXDIR}/mdrun -nt $nrpoc"
I want to modify the force field and topology files for SILCS simulation¶
As an example, if there is the need to add extra bonds that are not present in the standard force field definitions, this is the procedure to make the necessary modifications. For example, some proteins contain two metals ions adjacent to each other, in which case it may be useful to place a bond connecting the ions. The protein 3bi0 has two Zn ions adjacent to each other, and adding such a bond is useful to restrain the distance between the ions to that in the crystal structure. Please refer to the GROMACS documentation regarding to modify the .top and ffbonded.itp files.
First, run the following command. This will copy the basic force field
to and generate an initial topology file in 1_setup/, allowing
you to edit them.
$SILCSBIODIR/silcs/1_setup_silcs_boxes prot=<prot PDB>
Then edit the force field parameter file
1_setup/charmm36.ff/ffbonded.itp.
If you want to modify the topology file (e.g. to add an explicit bond
between the ions), copy 1_setup/<prot>_gmx.top.1.bak to
1_setup/<prot>_gmx.top. Then edit <prot>_gmx.top and add the
desired bond between the two ions in the [ bond ] list.
Once the files are edited, re-run the 1_setup command with the
skip_pdb2gmx=true keyword. This will preserve your edits and
create the necessary files to run the SILCS simulations.
$SILCSBIODIR/silcs/1_setup_silcs_boxes prot=<prot PDB> skip_pdb2gmx=true
Once this completes, run the $SILCSBIODIR/silcs/2a_run_gcmd script to
initiate your SILCS simulations.
I want to visualize FragMaps using MOE¶
By default, SILCS FragMaps are in the MAP grid file format. However, this file format is not supported in MOE. Please see FragMaps in MOE for detailed instructions on creating FragMaps in a MOE-compatible format.
How do I handle phosphorylated amino acids?¶
The following phosphorylated amino acids are supported:
pSer
pThr
pTyr
To create a phosphorylated amino acid, rename that amino acid in your input pdb file as follows:
SER => SP1 or SP2
THR => THP1 or THP2
TYR => TP1 or TP2
The number at the end of the amino acid name refers to whether the phosphate group has mono- or divalent charge.
What if my protein has a glycan attached to it?¶
While setting up a glycan-containing protein directly from a PDB file is not currently supported, you can set up your simulation system for SILCS if you have a PSF file created with the CHARMM36 force field.
An example can be found in the $SILCSBIODIR/examples/glycan folder.
Running the setup.sh script in that directory will run the example and
create a folder named 1_setup.
For your own system, copy the gromacs folder and setup.sh file
and edit the copied setup.sh file before running it:
psffile="psf/step1_pdbreader.psf" # PSF file
pdbfile="psf/step1_pdbreader.pdb" # PDB file
prefix="5vgp" # prefix for the SILCS simulation
What happens when I set up SILCS simulations with an input structure containing a metal ion?¶
SILCS simulation supports a variety of metal ions, including calcium, copper, iron, magnesium, manganese, nickel, and zinc. If an ion in the input structure is located close to protein atoms (~ 3A), the setup script will automatically create covalent linkages between the metal ion and nearby protein residues so as to ensure the coordination structure is maintained throughout the SILCS GCMC/MD simulations.
My protein contains iron and I want to set a +3 charge state¶
By default, the SILCS setup assigns a +2 charge to iron. If you want to change the charge of the iron ion, change the residue name of the ion in the input PDB to FE3. The setup script will then assign a +3 charge to that ion.
How do I fit my membrane protein in a bilayer as suggested by the OPM server?¶
Prepare your bare membrane protein with Alphafold, MOE etc. as <prot PDB>
Upload the <prot PDB> to the OPM server (PPM 3.0) and get the output.
Download the <OPM output pdb>, open it with PyMOL or VMD, and determine the required translation along the z-axis:
PyMOL
Select protein atoms:
sele all and polymerCalculate center of mass (COM):
centerofmass seleCopy the Z-coordinate_of_COM and go to the next step
VMD
Open VMD Main >> Extensions >> TkConsole
Create selection with protein atoms:
set sel1 [atomselect top “all and protein”]Calculate center of mass (COM):
set com1 [measure center $sel1 weight mass]Copy the Z-coordinate_of_COM and go to the next step
Align <prot PDB> to <OPM output pdb> for subsequent use OR extract protein from <OPM output pdb> for subsequent use.
Run setup with the additional
offset_z:$SILCSBIODIR/silcs-memb/1a_fit_protein_in_bilayer prot=<aligned prot PDB> orient_principal_axis=false offset_z=<Z-coordinate_of_COM>
How do I calculate the difference between two sets of SILCS FragMaps?¶
Copy the two silcs_fragmaps_xxx directories you wish to use to a new directory.
Please make sure:
Either both sets of FragMaps have the same grid setup, which means the header of the
silcs_fragmaps_xxx/maps/*.mapfiles are the same.Or both sets of FragMaps have were created relative to proteins having the same rotational and translational alignments. If the input proteins were not aligned to each other prior to running SILCS, you must use the “ref=” option with the
2b_gen_mapscommand to ensure alignment. A map cutting algorithm will be used to make the FragMaps compatible, with the smaller dimensions of the two grids used to calculate the difference maps.
Run the following command:
$SILCSBIODIR/utils/calc_difference_maps.sh [silcs_fragmaps_#1] [silcs_fragmaps_#2]
The difference maps will be saved in a directory
DIFF_MAPS_fragmap#1_VS_fragmap#2 and can be visualized in the same way
as regular FragMaps (Visualizing SILCS FragMaps).
Tip
In the difference maps, if the value at a grid point is negative, the grid free energy (GFE) of the grid point in fragmap_#1 is more favourable than in fragmap_#2 and if it is positive, then fragmap_#2 is more favourable than fragmap_#1.
Note
The difference maps are intended ONLY for visualization and NOT for any quantitative calculations.
How do I include a covalently bound ligand/cofactor in SILCS simulations?¶
While creating a covalent bond between your protein and the ligand/cofactor
during SILCS setup is not supported, it is possible to achieve a good
approximation of the structure by using positional restraints to maintain
relative geometries during the SILCS simulation. To do so, provide the
ligand/cofactor as an individual ligand using the option lig=lig.mol2
during setup. The mol2 file should contain only the ligand/cofactor molecule
and must have all hydrogen atoms and a appropriate three-dimensional
coordinates that provide a reasonable internal geometry and place it
correctly relative to the protein. Weak position restraints will be
automatically added on the non-hydrogen atoms of ligand/cofactor, so it will
be positionally restrained. (Note: you will also need to set the option
scramblesc=false when running 1_setup_silcs_boxes in order to
keep protein sideschain conformations from being scrambled.)
$SILCSBIODIR/silcs/1_setup_silcs_boxes prot=<prot PDB> lig=lig.mol2 scramblesc=false
However, the amino acid sidechain to which the ligand/cofactor is meant to be
covalently bound will be unrestrained. To address this, you will need to
manually add weak position restraints on non-hydrogen atoms of the sidechain.
To do this, edit the 1_setup/posre_protein_ca.itp file which restraints
C-alpha atoms and append it with the non-hydrogen atoms from your sidechain.
; position restraints for (atomname_CA_or_atomname_"C1'")
[ position_restraints ]
; i funct fcx fcy fcz
5 1 50.208 50.208 50.208
20 1 50.208 50.208 50.208
46 1 50.208 50.208 50.208
58 1 50.208 50.208 50.208
.
.
.
5571 1 50.208 50.208 50.208
5589 1 50.208 50.208 50.208
5603 1 50.208 50.208 50.208
5617 1 50.208 50.208 50.208
; **ADD YOUR SIDECHAIN ATOMS BELOW**
xxxx 1 50.208 50.208 50.208
xxxx 1 50.208 50.208 50.208
xxxx 1 50.208 50.208 50.208
xxxx 1 50.208 50.208 50.208
Please visualize the 1_setup/<prot>_silcs.1-10.pdb files
to make sure the ligand and sidechain are correctly prepared before you
run the 2a_run_gcmd command.