ISDB: setting up a Metadynamics Metainference simulation

Aims

The aim of this tutorial is to introduce the users to the ISDB module and in particular to Metadynamics Metainference [23] [24] ensemble determination. We will reproduce the setup of the simulation for a simple system [86] . For a general overview of the problem of ensembles determination please read [25] .

Objectives

Once this tutorial is completed students will be able to:

  • Setup their own PLUMED-ISDB simulation.

Resources

The TARBALL for this project contains the following files:

  • charmm36-eef1sb.ff: the force-field files for gromacs (not needed)
  • system: a folder with reference files for gromacs (not needed)
  • reference-impl: a folder to perform a simple implicit solvent simulation
  • reference-impl-pbmetad: a folder to perform a pbmetad implicit solvent simulation
  • m_and_m: a folder to perform a metadynamics metainference simulation

Introduction

Molecular dynamics simulations are the ideal tool to determine at atomistic resolution the behavior of complex molecules. This great resolution power comes at the cost of approximations that affects the agreement with actual experimental observables. At the same time experimental data alone are generally speaking not enough to determine a structural ensemble due the inverse nature of the problem, that is to go from few observables to many atoms in many different configurations. Furthermore, experimental data are affected by errors of multiple nature, from noise, systematic errors and errors in their atomistic interpretation. Most important experimental data are the result of the averaging over the ensemble of structure so it is not trivial to deconvolve this signal. One possibility is that of employing MD simulations together with experimental data to generate simulations already corrected for the data themselves. With METAINFERENCE this is done on-the-fly by adding an additional energy to the system that takes into account the agreement with the experimental data considering the multiple sources of errors.

Run a reference simulation

The system we use is the EGAAWAASS peptide used in ref. [86] . First of all we will run a simulation in implicit solvent using the EEF1-SB CHARMM36 force field. EEF1-SB includes a correction to the standard backbone torsion potential of CHARMM36, an electrostatic interaction with a distance dependent dielectric constant and a simple gaussian form for the solvation energy. The first two terms are implemented in the force field and using table potentials while the latter is implemented as a collective variable in PLUMED, EEFSOLV .

Click on the labels of the actions for more information on what each action computes
tested on v2.9
# this is optional and tell to VIM that this is a PLUMED file
# vim: ft=plumed
# see comments just below this input file
#SETTINGS MOLFILE=user-doc/tutorials/others/isdb-1/reference-impl/egaawaass.pdb
MOLINFO 
MOLTYPE
compulsory keyword ( default=protein )
what kind of molecule is contained in the pdb file - usually not needed since protein/RNA/DNA
are compatible 
=protein 
STRUCTURE
compulsory keyword 
a file in pdb format containing a reference structure. 
=egaawaass.pdb 
WHOLEMOLECULES 
ENTITY0
the atoms that make up a molecule that you wish to align. 
=1-111 
# EEF1SB Implicit solvation
#SETTINGS AUXFILE=user-doc/tutorials/others/isdb-1/reference-impl/index.ndx
protein-h: GROUP 
NDX_FILE
the name of index file (gromacs syntax) 
=index.ndx 
NDX_GROUP
the name of the group to be imported (gromacs syntax) - first group found is used
by default 
=Protein-H 
solv: EEFSOLV 
ATOMS
The atoms to be included in the calculation, e.g. 
=protein-h 
NOPBC
( default=off ) ignore the periodic boundary conditions when calculating distances

 
NL_STRIDE
compulsory keyword ( default=40 )
The frequency with which the neighbor list is updated. 
=20 
NL_BUFFER
compulsory keyword ( default=0.1 )
The buffer to the intrinsic cutoff used when calculating pairwise interactions. 
=0.1
bias: BIASVALUE ARG=solv

This can be run using gromacs (unfortunately recent versions of gromacs do not support Verlet groups with table potentials, so performances are currently sub-optimal on the gromacs side)

gmx_mpi mdrun -s run.tpr -table table.xvg -tablep table.xvg -plumed plumed-eef1.dat -v

In order to have a converged sampling for this reference ensemble calculation it is useful to setup a Metadynamics calculation. In particular we will use PBMETAD because it is then a natural choice for Metadynamics Metainference later.

Click on the labels of the actions for more information on what each action computes
tested on v2.9
# this is optional and tell to VIM that this is a PLUMED file
# vim: ft=plumed
# see comments just below this input file
#SETTINGS MOLFILE=user-doc/tutorials/others/isdb-1/reference-impl/egaawaass.pdb
MOLINFO 
MOLTYPE
compulsory keyword ( default=protein )
what kind of molecule is contained in the pdb file - usually not needed since protein/RNA/DNA
are compatible 
=protein 
STRUCTURE
compulsory keyword 
a file in pdb format containing a reference structure. 
=egaawaass.pdb 
WHOLEMOLECULES 
ENTITY0
the atoms that make up a molecule that you wish to align. 
=1-111 
# EEF1SB Implicit solvation
#SETTINGS AUXFILE=user-doc/tutorials/others/isdb-1/reference-impl/index.ndx
protein-h: GROUP 
NDX_FILE
the name of index file (gromacs syntax) 
=index.ndx 
NDX_GROUP
the name of the group to be imported (gromacs syntax) - first group found is used
by default 
=Protein-H 
solv: EEFSOLV 
ATOMS
The atoms to be included in the calculation, e.g. 
=protein-h 
NOPBC
( default=off ) ignore the periodic boundary conditions when calculating distances

 
NL_STRIDE
compulsory keyword ( default=40 )
The frequency with which the neighbor list is updated. 
=20 
NL_BUFFER
compulsory keyword ( default=0.1 )
The buffer to the intrinsic cutoff used when calculating pairwise interactions. 
=0.1 
bias: BIASVALUE 
ARG
the input for this action is the scalar output from one or more other actions. 
=solv 
# CVs, Psi9, Phi1 are not defined
psi1: TORSION 
ATOMS
the four atoms involved in the torsional angle 
=@psi-1 
NOPBC
( default=off ) ignore the periodic boundary conditions when calculating distances

 
psi2: TORSION 
ATOMS
the four atoms involved in the torsional angle 
=@psi-2 
NOPBC
( default=off ) ignore the periodic boundary conditions when calculating distances

 
psi3: TORSION 
ATOMS
the four atoms involved in the torsional angle 
=@psi-3 
NOPBC
( default=off ) ignore the periodic boundary conditions when calculating distances

 
psi4: TORSION 
ATOMS
the four atoms involved in the torsional angle 
=@psi-4 
NOPBC
( default=off ) ignore the periodic boundary conditions when calculating distances

 
psi5: TORSION 
ATOMS
the four atoms involved in the torsional angle 
=@psi-5 
NOPBC
( default=off ) ignore the periodic boundary conditions when calculating distances

 
psi6: TORSION 
ATOMS
the four atoms involved in the torsional angle 
=@psi-6 
NOPBC
( default=off ) ignore the periodic boundary conditions when calculating distances

 
psi7: TORSION 
ATOMS
the four atoms involved in the torsional angle 
=@psi-7 
NOPBC
( default=off ) ignore the periodic boundary conditions when calculating distances

 
psi8: TORSION 
ATOMS
the four atoms involved in the torsional angle 
=@psi-8 
NOPBC
( default=off ) ignore the periodic boundary conditions when calculating distances

 
phi2: TORSION 
ATOMS
the four atoms involved in the torsional angle 
=@phi-2 
NOPBC
( default=off ) ignore the periodic boundary conditions when calculating distances

 
phi3: TORSION 
ATOMS
the four atoms involved in the torsional angle 
=@phi-3 
NOPBC
( default=off ) ignore the periodic boundary conditions when calculating distances

 
phi4: TORSION 
ATOMS
the four atoms involved in the torsional angle 
=@phi-4 
NOPBC
( default=off ) ignore the periodic boundary conditions when calculating distances

 
phi5: TORSION 
ATOMS
the four atoms involved in the torsional angle 
=@phi-5 
NOPBC
( default=off ) ignore the periodic boundary conditions when calculating distances

 
phi6: TORSION 
ATOMS
the four atoms involved in the torsional angle 
=@phi-6 
NOPBC
( default=off ) ignore the periodic boundary conditions when calculating distances

 
phi7: TORSION 
ATOMS
the four atoms involved in the torsional angle 
=@phi-7 
NOPBC
( default=off ) ignore the periodic boundary conditions when calculating distances

 
phi8: TORSION 
ATOMS
the four atoms involved in the torsional angle 
=@phi-8 
NOPBC
( default=off ) ignore the periodic boundary conditions when calculating distances

 
phi9: TORSION 
ATOMS
the four atoms involved in the torsional angle 
=@phi-9 
NOPBC
( default=off ) ignore the periodic boundary conditions when calculating distances

 
ahc: ALPHARMSD 
RESIDUES
this command is used to specify the set of residues that could conceivably form part
of the secondary structure. 
=all 
TYPE
compulsory keyword ( default=DRMSD )
the manner in which RMSD alignment is performed. 
=OPTIMAL 
LESS_THAN
calculate the number of variables less than a certain target value. 
={RATIONAL R_0=0.12}  
# Bulky Trp residue dihedral
dihtrp_cacb: TORSION 
ATOMS
the four atoms involved in the torsional angle 
=67,47,49,52 
NOPBC
( default=off ) ignore the periodic boundary conditions when calculating distances

 
dihtrp_cbcg: TORSION 
ATOMS
the four atoms involved in the torsional angle 
=47,49,52,53 
NOPBC
( default=off ) ignore the periodic boundary conditions when calculating distances

 
protein-ca: GROUP 
NDX_FILE
the name of index file (gromacs syntax) 
=index.ndx 
NDX_GROUP
the name of the group to be imported (gromacs syntax) - first group found is used
by default 
=C-alpha 
gyr: GYRATION 
TYPE
compulsory keyword ( default=RADIUS )
The type of calculation relative to the Gyration Tensor you want to perform 
=RADIUS 
ATOMS
the group of atoms that you are calculating the Gyration Tensor for. 
=protein-ca 
NOPBC
( default=off ) ignore the periodic boundary conditions when calculating distances

 
# PBMetaD
pb: PBMETAD ...
   
   
ARG
the input for this action is the scalar output from one or more other actions. 
=phi2,phi3,phi4,phi5,phi6,phi7,phi8,phi9,psi1,psi2,psi3,psi4,psi5,psi6,psi7,psi8,dihtrp_cacb,dihtrp_cbcg,ahc.lessthan 
   
SIGMA
compulsory keyword 
the widths of the Gaussian hills 
=0.1,0.1,0.1,0.1,0.1,0.1,0.1,0.1,0.1,0.1,0.1,0.1,0.1,0.1,0.1,0.1,0.1,0.1,0.01 
   
HEIGHT
the height of the Gaussian hills, one for all biases. 
=0.5 
   
PACE
compulsory keyword 
the frequency for hill addition, one for all biases 
=400 
   
BIASFACTOR
use well tempered metadynamics with this bias factor, one for all biases. 
=20 
   
GRID_MIN
the lower bounds for the grid 
=-pi,-pi,-pi,-pi,-pi,-pi,-pi,-pi,-pi,-pi,-pi,-pi,-pi,-pi,-pi,-pi,-pi,-pi,0 
   
GRID_MAX
the upper bounds for the grid 
=pi,pi,pi,pi,pi,pi,pi,pi,pi,pi,pi,pi,pi,pi,pi,pi,pi,pi,5 
   
GRID_WSTRIDE
frequency for dumping the grid 
=5000 
   
WALKERS_MPI
( default=off ) Switch on MPI version of multiple walkers - not compatible with WALKERS_*
options other than WALKERS_DIR 
 
...
PRINT 
FILE
the name of the file on which to output these quantities 
=COLVAR 
ARG
the input for this action is the scalar output from one or more other actions. 
=phi2,phi3,phi4,phi5,phi6,phi7,phi8,phi9,psi1,psi2,psi3,psi4,psi5,psi6,psi7,psi8,dihtrp_cacb,dihtrp_cbcg,ahc.lessthan 
STRIDE
compulsory keyword ( default=1 )
the frequency with which the quantities of interest should be output 
=200 
PRINT 
FILE
the name of the file on which to output these quantities 
=ENERGY 
ARG
the input for this action is the scalar output from one or more other actions. 
=bias.bias,pb.bias 
STRIDE
compulsory keyword ( default=1 )
the frequency with which the quantities of interest should be output 
=200

In this case we are running a multiple-replica simulation where the sampling is used to parallelize the Metadynamics time-dependent potential through the use of multiple walkers.

mpiexec -np 14 gmx_mpi mdrun -s topolnew -multi 14 -plumed plumed-eef1-pbmetad.dat -table table.xvg -tablep table.xvg >& log.out &

Metadynamics Metainference

The former simulations should provide a converged (check for this) ensemble for the peptide. As shown in [86] the agreement with the multiple available NMR experimental data is not perfect. In order to generate an ensemble compatible with most of the available experimental data it is possible to include them in the simulation using METAINFERENCE . To do so the forward models for the data sets should be defined in the input file. In this case we have backbone chemical shifts, CS2BACKBONE ; residual dipolar couplings for two bonds, RDC ; and J-couplings for multiple atoms, JCOUPLING. Once the forward models are defined for the data sets, the calculated data together with the corresponding experimental values can be used to calculate the metainference score. The metainference score is additive so it can be split into multiple METAINFERENCE entries. In this case we are using two metainference entries for the two sets of RDCs because these are compared with the experimental data modulo a constant that should be unique each data set. Then we use one metainference for all the jcouplings and another one for the chemical shifts.

Click on the labels of the actions for more information on what each action computes
tested on v2.9
# this is optional and tell to VIM that this is a PLUMED file
# vim: ft=plumed
# see comments just below this input file
#SETTINGS NREPLICAS=2
#SETTINGS MOLFILE=user-doc/tutorials/others/isdb-1/reference-impl/egaawaass.pdb
MOLINFO 
MOLTYPE
compulsory keyword ( default=protein )
what kind of molecule is contained in the pdb file - usually not needed since protein/RNA/DNA
are compatible 
=protein 
STRUCTURE
compulsory keyword 
a file in pdb format containing a reference structure. 
=egaawaass.pdb 
WHOLEMOLECULES 
ENTITY0
the atoms that make up a molecule that you wish to align. 
=1-111 
# EEF1SB Implicit solvation
#SETTINGS AUXFILE=user-doc/tutorials/others/isdb-1/reference-impl/index.ndx
protein-h: GROUP 
NDX_FILE
the name of index file (gromacs syntax) 
=index.ndx 
NDX_GROUP
the name of the group to be imported (gromacs syntax) - first group found is used
by default 
=Protein-H 
solv: EEFSOLV 
ATOMS
The atoms to be included in the calculation, e.g. 
=protein-h 
NOPBC
( default=off ) ignore the periodic boundary conditions when calculating distances

 
NL_STRIDE
compulsory keyword ( default=40 )
The frequency with which the neighbor list is updated. 
=20 
NL_BUFFER
compulsory keyword ( default=0.1 )
The buffer to the intrinsic cutoff used when calculating pairwise interactions. 
=0.1 
bias: BIASVALUE 
ARG
the input for this action is the scalar output from one or more other actions. 
=solv 
# CVs, Psi9, Phi1 are not defined
psi1: TORSION 
ATOMS
the four atoms involved in the torsional angle 
=@psi-1 
NOPBC
( default=off ) ignore the periodic boundary conditions when calculating distances

 
psi2: TORSION 
ATOMS
the four atoms involved in the torsional angle 
=@psi-2 
NOPBC
( default=off ) ignore the periodic boundary conditions when calculating distances

 
psi3: TORSION 
ATOMS
the four atoms involved in the torsional angle 
=@psi-3 
NOPBC
( default=off ) ignore the periodic boundary conditions when calculating distances

 
psi4: TORSION 
ATOMS
the four atoms involved in the torsional angle 
=@psi-4 
NOPBC
( default=off ) ignore the periodic boundary conditions when calculating distances

 
psi5: TORSION 
ATOMS
the four atoms involved in the torsional angle 
=@psi-5 
NOPBC
( default=off ) ignore the periodic boundary conditions when calculating distances

 
psi6: TORSION 
ATOMS
the four atoms involved in the torsional angle 
=@psi-6 
NOPBC
( default=off ) ignore the periodic boundary conditions when calculating distances

 
psi7: TORSION 
ATOMS
the four atoms involved in the torsional angle 
=@psi-7 
NOPBC
( default=off ) ignore the periodic boundary conditions when calculating distances

 
psi8: TORSION 
ATOMS
the four atoms involved in the torsional angle 
=@psi-8 
NOPBC
( default=off ) ignore the periodic boundary conditions when calculating distances

 
phi2: TORSION 
ATOMS
the four atoms involved in the torsional angle 
=@phi-2 
NOPBC
( default=off ) ignore the periodic boundary conditions when calculating distances

 
phi3: TORSION 
ATOMS
the four atoms involved in the torsional angle 
=@phi-3 
NOPBC
( default=off ) ignore the periodic boundary conditions when calculating distances

 
phi4: TORSION 
ATOMS
the four atoms involved in the torsional angle 
=@phi-4 
NOPBC
( default=off ) ignore the periodic boundary conditions when calculating distances

 
phi5: TORSION 
ATOMS
the four atoms involved in the torsional angle 
=@phi-5 
NOPBC
( default=off ) ignore the periodic boundary conditions when calculating distances

 
phi6: TORSION 
ATOMS
the four atoms involved in the torsional angle 
=@phi-6 
NOPBC
( default=off ) ignore the periodic boundary conditions when calculating distances

 
phi7: TORSION 
ATOMS
the four atoms involved in the torsional angle 
=@phi-7 
NOPBC
( default=off ) ignore the periodic boundary conditions when calculating distances

 
phi8: TORSION 
ATOMS
the four atoms involved in the torsional angle 
=@phi-8 
NOPBC
( default=off ) ignore the periodic boundary conditions when calculating distances

 
phi9: TORSION 
ATOMS
the four atoms involved in the torsional angle 
=@phi-9 
NOPBC
( default=off ) ignore the periodic boundary conditions when calculating distances

 
ahc: ALPHARMSD 
RESIDUES
this command is used to specify the set of residues that could conceivably form part
of the secondary structure. 
=all 
TYPE
compulsory keyword ( default=DRMSD )
the manner in which RMSD alignment is performed. 
=OPTIMAL 
LESS_THAN
calculate the number of variables less than a certain target value. 
={RATIONAL R_0=0.12}  
# Bulky Trp residue dihedral
dihtrp_cacb: TORSION 
ATOMS
the four atoms involved in the torsional angle 
=67,47,49,52 
NOPBC
( default=off ) ignore the periodic boundary conditions when calculating distances

 
dihtrp_cbcg: TORSION 
ATOMS
the four atoms involved in the torsional angle 
=47,49,52,53 
NOPBC
( default=off ) ignore the periodic boundary conditions when calculating distances

 
protein-ca: GROUP 
NDX_FILE
the name of index file (gromacs syntax) 
=index.ndx 
NDX_GROUP
the name of the group to be imported (gromacs syntax) - first group found is used
by default 
=C-alpha 
gyr: GYRATION 
TYPE
compulsory keyword ( default=RADIUS )
The type of calculation relative to the Gyration Tensor you want to perform 
=RADIUS 
ATOMS
the group of atoms that you are calculating the Gyration Tensor for. 
=protein-ca 
NOPBC
( default=off ) ignore the periodic boundary conditions when calculating distances

 
# PBMetaD
pb: PBMETAD ...
   
   
ARG
the input for this action is the scalar output from one or more other actions. 
=phi2,phi3,phi4,phi5,phi6,phi7,phi8,phi9,psi1,psi2,psi3,psi4,psi5,psi6,psi7,psi8,dihtrp_cacb,dihtrp_cbcg,ahc.lessthan 
   
SIGMA
compulsory keyword 
the widths of the Gaussian hills 
=0.1,0.1,0.1,0.1,0.1,0.1,0.1,0.1,0.1,0.1,0.1,0.1,0.1,0.1,0.1,0.1,0.1,0.1,0.01 
   
HEIGHT
the height of the Gaussian hills, one for all biases. 
=0.5 
   
PACE
compulsory keyword 
the frequency for hill addition, one for all biases 
=400 
   
BIASFACTOR
use well tempered metadynamics with this bias factor, one for all biases. 
=20 
   
GRID_MIN
the lower bounds for the grid 
=-pi,-pi,-pi,-pi,-pi,-pi,-pi,-pi,-pi,-pi,-pi,-pi,-pi,-pi,-pi,-pi,-pi,-pi,0 
   
GRID_MAX
the upper bounds for the grid 
=pi,pi,pi,pi,pi,pi,pi,pi,pi,pi,pi,pi,pi,pi,pi,pi,pi,pi,5 
   
GRID_WSTRIDE
frequency for dumping the grid 
=5000 
   
WALKERS_MPI
( default=off ) Switch on MPI version of multiple walkers - not compatible with WALKERS_*
options other than WALKERS_DIR 
 
...
PRINT 
FILE
the name of the file on which to output these quantities 
=COLVAR 
ARG
the input for this action is the scalar output from one or more other actions. 
=phi2,phi3,phi4,phi5,phi6,phi7,phi8,phi9,psi1,psi2,psi3,psi4,psi5,psi6,psi7,psi8,dihtrp_cacb,dihtrp_cbcg,ahc.lessthan 
STRIDE
compulsory keyword ( default=1 )
the frequency with which the quantities of interest should be output 
=200 
PRINT 
FILE
the name of the file on which to output these quantities 
=ENERGY 
ARG
the input for this action is the scalar output from one or more other actions. 
=bias.bias,pb.bias 
STRIDE
compulsory keyword ( default=1 )
the frequency with which the quantities of interest should be output 
=200 
# EXPERIMENTAL DATA SECTION
# RDCs (Grzesiek et al.)
# xGAAWAASS
nh: RDC ...
   
GYROM
compulsory keyword ( default=1. )
Add the product of the gyromagnetic constants for the bond. 
=-72.5388 
   
SCALE
compulsory keyword ( default=1. )
Add the scaling factor to take into account concentration and other effects. 
=0.0001 
   
NOPBC
( default=off ) ignore the periodic boundary conditions when calculating distances

 
   
ATOMS1
the couple of atoms involved in each of the bonds for which you wish to calculate
the RDC. 
=18,19 
COUPLING1
Add an experimental value for each coupling (needed by SVD and useful for \ref STATS)..

=-5.4 
   
ATOMS2
the couple of atoms involved in each of the bonds for which you wish to calculate
the RDC. 
=25,26 
COUPLING2
Add an experimental value for each coupling (needed by SVD and useful for \ref STATS)..

=-1.26 
   
ATOMS3
the couple of atoms involved in each of the bonds for which you wish to calculate
the RDC. 
=35,36 
COUPLING3
Add an experimental value for each coupling (needed by SVD and useful for \ref STATS)..

=-5.22 
   
ATOMS4
the couple of atoms involved in each of the bonds for which you wish to calculate
the RDC. 
=45,46 
COUPLING4
Add an experimental value for each coupling (needed by SVD and useful for \ref STATS)..

=-0.91 
   
ATOMS5
the couple of atoms involved in each of the bonds for which you wish to calculate
the RDC. 
=69,70 
COUPLING5
Add an experimental value for each coupling (needed by SVD and useful for \ref STATS)..

=2.33 
   
ATOMS6
the couple of atoms involved in each of the bonds for which you wish to calculate
the RDC. 
=79,80 
COUPLING6
Add an experimental value for each coupling (needed by SVD and useful for \ref STATS)..

=-2.88 
   
ATOMS7
the couple of atoms involved in each of the bonds for which you wish to calculate
the RDC. 
=89,90 
COUPLING7
Add an experimental value for each coupling (needed by SVD and useful for \ref STATS)..

=-8.37 
   
ATOMS8
the couple of atoms involved in each of the bonds for which you wish to calculate
the RDC. 
=100,101 
COUPLING8
Add an experimental value for each coupling (needed by SVD and useful for \ref STATS)..

=-3.78 
   
...
# ExAAWAASx
caha: RDC ...
   
GYROM
compulsory keyword ( default=1. )
Add the product of the gyromagnetic constants for the bond. 
=179.9319 
   
SCALE
compulsory keyword ( default=1. )
Add the scaling factor to take into account concentration and other effects. 
=0.0001 
   
NOPBC
( default=off ) ignore the periodic boundary conditions when calculating distances

 
   
ATOMS1
the couple of atoms involved in each of the bonds for which you wish to calculate
the RDC. 
=5,6 
COUPLING1
Add an experimental value for each coupling (needed by SVD and useful for \ref STATS)..

=12.95 
   
ATOMS2
the couple of atoms involved in each of the bonds for which you wish to calculate
the RDC. 
=27,28 
COUPLING2
Add an experimental value for each coupling (needed by SVD and useful for \ref STATS)..

=11.5 
   
ATOMS3
the couple of atoms involved in each of the bonds for which you wish to calculate
the RDC. 
=37,38 
COUPLING3
Add an experimental value for each coupling (needed by SVD and useful for \ref STATS)..

=21.42 
   
ATOMS4
the couple of atoms involved in each of the bonds for which you wish to calculate
the RDC. 
=47,48 
COUPLING4
Add an experimental value for each coupling (needed by SVD and useful for \ref STATS)..

=-9.37 
   
ATOMS5
the couple of atoms involved in each of the bonds for which you wish to calculate
the RDC. 
=71,72 
COUPLING5
Add an experimental value for each coupling (needed by SVD and useful for \ref STATS)..

=10.01 
   
ATOMS6
the couple of atoms involved in each of the bonds for which you wish to calculate
the RDC. 
=81,82 
COUPLING6
Add an experimental value for each coupling (needed by SVD and useful for \ref STATS)..

=15.01 
   
ATOMS7
the couple of atoms involved in each of the bonds for which you wish to calculate
the RDC. 
=91,92 
COUPLING7
Add an experimental value for each coupling (needed by SVD and useful for \ref STATS)..

=15.73 
   
...
#RDCS
byrdcnh: METAINFERENCE ...
   
ARG
the input for this action is the scalar output from one or more other actions. 
=(nh\.rdc-.*),pb.bias 
   
PARARG
reference values for the experimental data, these can be provided as arguments without
derivatives 
=(nh\.exp-.*) 
   
REWEIGHT
( default=off ) simple REWEIGHT using the latest ARG as energy 
 
   
NOISETYPE
compulsory keyword ( default=MGAUSS )
functional form of the noise (GAUSS,MGAUSS,OUTLIERS,MOUTLIERS,GENERIC) 
=MGAUSS 
   
OPTSIGMAMEAN
compulsory keyword ( default=NONE )
Set to NONE/SEM to manually set sigma mean, or to estimate it on the fly 
=SEM_MAX 
AVERAGING
Stride for calculation of averaged weights and sigma_mean 
=400 
   
SCALEDATA
( default=off ) Set to TRUE if you want to sample a scaling factor common to all
values and replicas 
 
SCALE_PRIOR
compulsory keyword ( default=FLAT )
either FLAT or GAUSSIAN 
=GAUSSIAN 
SCALE0
 could not find this keyword 
=8.0 
DSCALE
maximum MC move of the scaling factor 
=0.5 
   
SIGMA0
 could not find this keyword 
=5.0 
SIGMA_MIN
compulsory keyword ( default=0.0 )
minimum value of the uncertainty parameter 
=0.0001 
SIGMA_MAX
compulsory keyword ( default=10. )
maximum value of the uncertainty parameter 
=15.0 
DSIGMA
maximum MC move of the uncertainty parameter 
=0.1 
   
WRITE_STRIDE
compulsory keyword ( default=10000 )
write the status to a file every N steps, this can be used for restart/continuation

=10000 
   
...
#RDCS
byrdccaha: METAINFERENCE ...
   
ARG
the input for this action is the scalar output from one or more other actions. 
=(caha\.rdc-.*),pb.bias 
   
PARARG
reference values for the experimental data, these can be provided as arguments without
derivatives 
=(caha\.exp-.*) 
   
REWEIGHT
( default=off ) simple REWEIGHT using the latest ARG as energy 
 
   
NOISETYPE
compulsory keyword ( default=MGAUSS )
functional form of the noise (GAUSS,MGAUSS,OUTLIERS,MOUTLIERS,GENERIC) 
=MGAUSS 
   
OPTSIGMAMEAN
compulsory keyword ( default=NONE )
Set to NONE/SEM to manually set sigma mean, or to estimate it on the fly 
=SEM_MAX 
AVERAGING
Stride for calculation of averaged weights and sigma_mean 
=400 
   
SCALEDATA
( default=off ) Set to TRUE if you want to sample a scaling factor common to all
values and replicas 
 
SCALE_PRIOR
compulsory keyword ( default=FLAT )
either FLAT or GAUSSIAN 
=GAUSSIAN 
SCALE0
 could not find this keyword 
=9.0 
DSCALE
maximum MC move of the scaling factor 
=0.5 
   
SIGMA0
 could not find this keyword 
=5.0 
SIGMA_MIN
compulsory keyword ( default=0.0 )
minimum value of the uncertainty parameter 
=0.0001 
SIGMA_MAX
compulsory keyword ( default=10. )
maximum value of the uncertainty parameter 
=15.0 
DSIGMA
maximum MC move of the uncertainty parameter 
=0.1 
   
WRITE_STRIDE
compulsory keyword ( default=10000 )
write the status to a file every N steps, this can be used for restart/continuation

=10000 
   
...
# xGxAWxASx
jhan: JCOUPLING ...
   
TYPE
compulsory keyword 
Type of J-coupling to compute (HAN,HAHN,CCG,NCG,CUSTOM) 
=HAN 
   
NOPBC
( default=off ) ignore the periodic boundary conditions when calculating distances

 
   
ATOMS1
the 4 atoms involved in each of the bonds for which you wish to calculate the J-coupling.

=@psi-2 
COUPLING1
Add an experimental value for each coupling. 
=-0.49 
   
ATOMS2
the 4 atoms involved in each of the bonds for which you wish to calculate the J-coupling.

=@psi-4 
COUPLING2
Add an experimental value for each coupling. 
=-0.54 
   
ATOMS3
the 4 atoms involved in each of the bonds for which you wish to calculate the J-coupling.

=@psi-5 
COUPLING3
Add an experimental value for each coupling. 
=-0.53 
   
ATOMS4
the 4 atoms involved in each of the bonds for which you wish to calculate the J-coupling.

=@psi-7 
COUPLING4
Add an experimental value for each coupling. 
=-0.39 
   
ATOMS5
the 4 atoms involved in each of the bonds for which you wish to calculate the J-coupling.

=@psi-8 
COUPLING5
Add an experimental value for each coupling. 
=-0.39 
   
...
# xxAAWAASS
jhahn: JCOUPLING ...
   
TYPE
compulsory keyword 
Type of J-coupling to compute (HAN,HAHN,CCG,NCG,CUSTOM) 
=HAHN 
   
NOPBC
( default=off ) ignore the periodic boundary conditions when calculating distances

 
   
ATOMS1
the 4 atoms involved in each of the bonds for which you wish to calculate the J-coupling.

=@phi-2 
COUPLING1
Add an experimental value for each coupling. 
=6.05 
   
ATOMS2
the 4 atoms involved in each of the bonds for which you wish to calculate the J-coupling.

=@phi-3 
COUPLING2
Add an experimental value for each coupling. 
=5.95 
   
ATOMS3
the 4 atoms involved in each of the bonds for which you wish to calculate the J-coupling.

=@phi-4 
COUPLING3
Add an experimental value for each coupling. 
=6.44 
   
ATOMS4
the 4 atoms involved in each of the bonds for which you wish to calculate the J-coupling.

=@phi-5 
COUPLING4
Add an experimental value for each coupling. 
=6.53 
   
ATOMS5
the 4 atoms involved in each of the bonds for which you wish to calculate the J-coupling.

=@phi-6 
COUPLING5
Add an experimental value for each coupling. 
=5.93 
   
ATOMS6
the 4 atoms involved in each of the bonds for which you wish to calculate the J-coupling.

=@phi-7 
COUPLING6
Add an experimental value for each coupling. 
=6.98 
   
ATOMS7
the 4 atoms involved in each of the bonds for which you wish to calculate the J-coupling.

=@phi-8 
COUPLING7
Add an experimental value for each coupling. 
=7.16 
   
...
# xxxxWxxxx
jccg: JCOUPLING ...
   
TYPE
compulsory keyword 
Type of J-coupling to compute (HAN,HAHN,CCG,NCG,CUSTOM) 
=CCG 
   
NOPBC
( default=off ) ignore the periodic boundary conditions when calculating distances

 
   
ATOMS1
the 4 atoms involved in each of the bonds for which you wish to calculate the J-coupling.

=@chi1-5 
COUPLING1
Add an experimental value for each coupling. 
=1.59 
   
...
# xxxxWxxxx
jncg: JCOUPLING ...
   
TYPE
compulsory keyword 
Type of J-coupling to compute (HAN,HAHN,CCG,NCG,CUSTOM) 
=NCG 
   
NOPBC
( default=off ) ignore the periodic boundary conditions when calculating distances

 
   
ATOMS1
the 4 atoms involved in each of the bonds for which you wish to calculate the J-coupling.

=@chi1-5 
COUPLING1
Add an experimental value for each coupling. 
=1.21 
   
...
#JC METAINFERENCE
byj: METAINFERENCE ...
   
ARG
the input for this action is the scalar output from one or more other actions. 
=(jhan\.j-.*),(jhahn\.j-.*),(jccg\.j.*),(jncg\.j.*),pb.bias 
   
PARARG
reference values for the experimental data, these can be provided as arguments without
derivatives 
=(jhan\.exp-.*),(jhahn\.exp-.*),(jccg\.exp.*),(jncg\.exp.*) 
   
REWEIGHT
( default=off ) simple REWEIGHT using the latest ARG as energy 
 
   
NOISETYPE
compulsory keyword ( default=MGAUSS )
functional form of the noise (GAUSS,MGAUSS,OUTLIERS,MOUTLIERS,GENERIC) 
=MGAUSS 
   
OPTSIGMAMEAN
compulsory keyword ( default=NONE )
Set to NONE/SEM to manually set sigma mean, or to estimate it on the fly 
=SEM_MAX 
AVERAGING
Stride for calculation of averaged weights and sigma_mean 
=400 
   
SIGMA0
 could not find this keyword 
=5.0 
SIGMA_MIN
compulsory keyword ( default=0.0 )
minimum value of the uncertainty parameter 
=0.0001 
SIGMA_MAX
compulsory keyword ( default=10. )
maximum value of the uncertainty parameter 
=15.0 
DSIGMA
maximum MC move of the uncertainty parameter 
=0.1 
   
WRITE_STRIDE
compulsory keyword ( default=10000 )
write the status to a file every N steps, this can be used for restart/continuation

=10000 
   
...

# Chemical shifts
#SETTINGS AUXFOLDER=user-doc/tutorials/others/isdb-1/m_and_m/data
# Chemical shifts
cs: CS2BACKBONE ...
   
NOPBC
( default=off ) ignore the periodic boundary conditions when calculating distances

 
ATOMS
The atoms to be included in the calculation, e.g. 
=1-111 
DATADIR
compulsory keyword ( default=data/ )
The folder with the experimental chemical shifts. 
=data 
TEMPLATE
compulsory keyword ( default=template.pdb )
A PDB file of the protein system. 
=egaawaass.pdb 
DOSCORE
( default=off ) activate metainference 
 
ARG
the input for this action is the scalar output from one or more other actions. 
=pb.bias 
   
NOISETYPE
compulsory keyword ( default=MGAUSS )
functional form of the noise (GAUSS,MGAUSS,OUTLIERS,MOUTLIERS,GENERIC) 
=MGAUSS 
REWEIGHT
( default=off ) simple REWEIGHT using the ARG as energy 
 
OPTSIGMAMEAN
compulsory keyword ( default=NONE )
Set to NONE/SEM to manually set sigma mean, or to estimate it on the fly 
=SEM_MAX 
AVERAGING
Stride for calculation of averaged weights and sigma_mean 
=400 
   
SIGMA0
 could not find this keyword 
=4.0 
   
SIGMA_MIN
compulsory keyword ( default=0.0 )
minimum value of the uncertainty parameter 
=0.0001 
   
SIGMA_MAX
compulsory keyword ( default=10. )
maximum value of the uncertainty parameter 
=5.00 
   
WRITE_STRIDE
compulsory keyword ( default=10000 )
write the status to a file every N steps, this can be used for restart/continuation

=10000 
...
mcs: BIASVALUE 
ARG
the input for this action is the scalar output from one or more other actions. 
=cs.score 
# output from METAINFERENCE
PRINT 
ARG
the input for this action is the scalar output from one or more other actions. 
=byrdcnh.* 
STRIDE
compulsory keyword ( default=1 )
the frequency with which the quantities of interest should be output 
=200 
FILE
the name of the file on which to output these quantities 
=BAYES.RDC.NH 
PRINT 
ARG
the input for this action is the scalar output from one or more other actions. 
=byrdccaha.* 
STRIDE
compulsory keyword ( default=1 )
the frequency with which the quantities of interest should be output 
=200 
FILE
the name of the file on which to output these quantities 
=BAYES.RDC.CAHA 
PRINT 
ARG
the input for this action is the scalar output from one or more other actions. 
=byj.* 
STRIDE
compulsory keyword ( default=1 )
the frequency with which the quantities of interest should be output 
=200 
FILE
the name of the file on which to output these quantities 
=BAYES.J 
PRINT 
ARG
the input for this action is the scalar output from one or more other actions. 
=cs.* 
STRIDE
compulsory keyword ( default=1 )
the frequency with which the quantities of interest should be output 
=200 
FILE
the name of the file on which to output these quantities 
=BAYES.CS 
# the following are useful for the analysis on-the-fly of the quality of the agreement with the experimentl data
ens: ENSEMBLE ...
   
ARG
the input for this action is the scalar output from one or more other actions. 
=(nh\.rdc-.*),(caha\.rdc-.*),(jhan\.j-.*),(jhahn\.j-.*),(jccg\.j-.*),(jncg\.j-.*),(cs\...-.*),pb.bias 
REWEIGHT
( default=off ) simple REWEIGHT using the latest ARG as energy 
 
   
...
nhst: STATS ...
   
ARG
the input for this action is the scalar output from one or more other actions. 
=(ens\.nh\.rdc-.*) 
PARARG
the input for this action is the scalar output from one or more other actions without
derivatives. 
=(nh\.exp-.*) 
   
...
cahast: STATS ...
   
ARG
the input for this action is the scalar output from one or more other actions. 
=(ens\.caha\.rdc-.*) 
PARARG
the input for this action is the scalar output from one or more other actions without
derivatives. 
=(caha\.exp-.*) 
   
...
csst: STATS ...
   
ARG
the input for this action is the scalar output from one or more other actions. 
=(ens\.cs\...-.*) 
PARARG
the input for this action is the scalar output from one or more other actions without
derivatives. 
=(cs\.exp.*) 
   
...
jhanst: STATS ...
   
ARG
the input for this action is the scalar output from one or more other actions. 
=(ens\.jhan\.j-.*) 
PARARG
the input for this action is the scalar output from one or more other actions without
derivatives. 
=(jhan\.exp-.*) 
   
...
jhahnst: STATS ...
   
ARG
the input for this action is the scalar output from one or more other actions. 
=(ens\.jhahn\.j-.*) 
PARARG
the input for this action is the scalar output from one or more other actions without
derivatives. 
=(jhahn\.exp-.*) 
   
...
jw5ccyst: STATS ...
   
ARG
the input for this action is the scalar output from one or more other actions. 
=(ens\.jccg\.j.*),(ens\.jccg\.j.*) 
PARARG
the input for this action is the scalar output from one or more other actions without
derivatives. 
=(jccg\.exp-.*),(jccg\.exp-.*) 
   
SQDEVSUM
( default=off ) calculates only SQDEVSUM 
 
   
...
jw5ncyst: STATS ...
   
ARG
the input for this action is the scalar output from one or more other actions. 
=(ens\.jncg\.j.*),(ens\.jncg\.j.*) 
PARARG
the input for this action is the scalar output from one or more other actions without
derivatives. 
=(jncg\.exp-.*),(jncg\.exp-.*) 
   
SQDEVSUM
( default=off ) calculates only SQDEVSUM 
 
   
...
#output from STATS
PRINT 
ARG
the input for this action is the scalar output from one or more other actions. 
=nhst.* 
STRIDE
compulsory keyword ( default=1 )
the frequency with which the quantities of interest should be output 
=2000 
FILE
the name of the file on which to output these quantities 
=ST.RDC.NH 
PRINT 
ARG
the input for this action is the scalar output from one or more other actions. 
=cahast.* 
STRIDE
compulsory keyword ( default=1 )
the frequency with which the quantities of interest should be output 
=2000 
FILE
the name of the file on which to output these quantities 
=ST.RDC.CAHA 
PRINT 
ARG
the input for this action is the scalar output from one or more other actions. 
=csst.* 
STRIDE
compulsory keyword ( default=1 )
the frequency with which the quantities of interest should be output 
=2000 
FILE
the name of the file on which to output these quantities 
=ST.CS
PRINT ARG=jhanst.*,jhahnst.*,jw5ccyst.*,jw5ncyst.* STRIDE=2000 FILE=ST.J

As for the former case we are running a multiple-replica simulation where in addition to multiple-walker metadynamics we are also coupling the replicas through Metainference. The use of multiple-walkers metadynamics is here key in order to have the same bias defined for all the replicas. This allows us to calculate a weighted average of the experimental observables where the weights are defined univocally from the bias [24] .

mpiexec -np 14 gmx_mpi mdrun -s topolnew -multi 14 -plumed plumed-eef1-pbmetad-m_m.dat -table table.xvg -tablep table.xvg >& log.out &