The aim of this tutorial is to introduce the users to the ISDB module and in particular to Metadynamics Metainference [18] [19] ensemble determination. We will reproduce the setup of the simulation for a simple system [73] . For a general overview of the problem of ensembles determination please read [20] .
Once this tutorial is completed students will be able to:
The TARBALL for this project contains the following files:
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.
The system we use is the EGAAWAASS peptide used in ref. [73] . 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 .
# 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 MOLINFObias: BIASVALUE ARG=solvMOLTYPE=proteincompulsory keyword ( default=protein ) what kind of molecule is contained in the pdb file - usually not needed since protein/RNA/DNA are compatibleSTRUCTURE=egaawaass.pdb WHOLEMOLECULEScompulsory keyword a file in pdb format containing a reference structure.ENTITY0=1-111 # EEF1SB Implicit solvation #SETTINGS AUXFILE=user-doc/tutorials/others/isdb-1/reference-impl/index.ndx protein-h: GROUPthe atoms that make up a molecule that you wish to align.NDX_FILE=index.ndxthe name of index file (gromacs syntax)NDX_GROUP=Protein-H solv: EEFSOLVthe name of the group to be imported (gromacs syntax) - first group found is used by defaultATOMS=protein-hThe atoms to be included in the calculation, e.g.NOPBC( default=off ) ignore the periodic boundary conditions when calculating distancesNL_STRIDE=20compulsory keyword ( default=40 ) The frequency with which the neighbor list is updated.NL_BUFFER=0.1compulsory keyword ( default=0.1 ) The buffer to the intrinsic cutoff used when calculating pairwise interactions.
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.
# 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 MOLINFOMOLTYPE=proteincompulsory keyword ( default=protein ) what kind of molecule is contained in the pdb file - usually not needed since protein/RNA/DNA are compatibleSTRUCTURE=egaawaass.pdb WHOLEMOLECULEScompulsory keyword a file in pdb format containing a reference structure.ENTITY0=1-111 # EEF1SB Implicit solvation #SETTINGS AUXFILE=user-doc/tutorials/others/isdb-1/reference-impl/index.ndx protein-h: GROUPthe atoms that make up a molecule that you wish to align.NDX_FILE=index.ndxthe name of index file (gromacs syntax)NDX_GROUP=Protein-H solv: EEFSOLVthe name of the group to be imported (gromacs syntax) - first group found is used by defaultATOMS=protein-hThe atoms to be included in the calculation, e.g.NOPBC( default=off ) ignore the periodic boundary conditions when calculating distancesNL_STRIDE=20compulsory keyword ( default=40 ) The frequency with which the neighbor list is updated.NL_BUFFER=0.1 bias: BIASVALUEcompulsory keyword ( default=0.1 ) The buffer to the intrinsic cutoff used when calculating pairwise interactions.ARG=solv # CVs, Psi9, Phi1 are not defined psi1: TORSIONthe input for this action is the scalar output from one or more other actions.ATOMS=@psi-1the four atoms involved in the torsional angleNOPBCpsi2: TORSION( default=off ) ignore the periodic boundary conditions when calculating distancesATOMS=@psi-2the four atoms involved in the torsional angleNOPBCpsi3: TORSION( default=off ) ignore the periodic boundary conditions when calculating distancesATOMS=@psi-3the four atoms involved in the torsional angleNOPBCpsi4: TORSION( default=off ) ignore the periodic boundary conditions when calculating distancesATOMS=@psi-4the four atoms involved in the torsional angleNOPBCpsi5: TORSION( default=off ) ignore the periodic boundary conditions when calculating distancesATOMS=@psi-5the four atoms involved in the torsional angleNOPBCpsi6: TORSION( default=off ) ignore the periodic boundary conditions when calculating distancesATOMS=@psi-6the four atoms involved in the torsional angleNOPBCpsi7: TORSION( default=off ) ignore the periodic boundary conditions when calculating distancesATOMS=@psi-7the four atoms involved in the torsional angleNOPBCpsi8: TORSION( default=off ) ignore the periodic boundary conditions when calculating distancesATOMS=@psi-8the four atoms involved in the torsional angleNOPBCphi2: TORSION( default=off ) ignore the periodic boundary conditions when calculating distancesATOMS=@phi-2the four atoms involved in the torsional angleNOPBCphi3: TORSION( default=off ) ignore the periodic boundary conditions when calculating distancesATOMS=@phi-3the four atoms involved in the torsional angleNOPBCphi4: TORSION( default=off ) ignore the periodic boundary conditions when calculating distancesATOMS=@phi-4the four atoms involved in the torsional angleNOPBCphi5: TORSION( default=off ) ignore the periodic boundary conditions when calculating distancesATOMS=@phi-5the four atoms involved in the torsional angleNOPBCphi6: TORSION( default=off ) ignore the periodic boundary conditions when calculating distancesATOMS=@phi-6the four atoms involved in the torsional angleNOPBCphi7: TORSION( default=off ) ignore the periodic boundary conditions when calculating distancesATOMS=@phi-7the four atoms involved in the torsional angleNOPBCphi8: TORSION( default=off ) ignore the periodic boundary conditions when calculating distancesATOMS=@phi-8the four atoms involved in the torsional angleNOPBCphi9: TORSION( default=off ) ignore the periodic boundary conditions when calculating distancesATOMS=@phi-9the four atoms involved in the torsional angleNOPBCahc: ALPHARMSD( default=off ) ignore the periodic boundary conditions when calculating distancesRESIDUES=allthis command is used to specify the set of residues that could conceivably form part of the secondary structure.TYPE=OPTIMALcompulsory keyword ( default=DRMSD ) the manner in which RMSD alignment is performed.LESS_THAN={RATIONAL R_0=0.12} # Bulky Trp residue dihedral dihtrp_cacb: TORSIONcalculate the number of variables less than a certain target value.ATOMS=67,47,49,52the four atoms involved in the torsional angleNOPBCdihtrp_cbcg: TORSION( default=off ) ignore the periodic boundary conditions when calculating distancesATOMS=47,49,52,53the four atoms involved in the torsional angleNOPBCprotein-ca: GROUP( default=off ) ignore the periodic boundary conditions when calculating distancesNDX_FILE=index.ndxthe name of index file (gromacs syntax)NDX_GROUP=C-alpha gyr: GYRATIONthe name of the group to be imported (gromacs syntax) - first group found is used by defaultTYPE=RADIUScompulsory keyword ( default=RADIUS ) The type of calculation relative to the Gyration Tensor you want to performATOMS=protein-cathe group of atoms that you are calculating the Gyration Tensor for.NOPBC# PBMetaD pb: PBMETAD ...( default=off ) ignore the periodic boundary conditions when calculating distancesARG=phi2,phi3,phi4,phi5,phi6,phi7,phi8,phi9,psi1,psi2,psi3,psi4,psi5,psi6,psi7,psi8,dihtrp_cacb,dihtrp_cbcg,ahc.lessthanthe input for this action is the scalar output from one or more other actions.SIGMA=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.01compulsory keyword the widths of the Gaussian hillsHEIGHT=0.5the height of the Gaussian hills, one for all biases.PACE=400compulsory keyword the frequency for hill addition, one for all biasesBIASFACTOR=20use well tempered metadynamics with this bias factor, one for all biases.GRID_MIN=-pi,-pi,-pi,-pi,-pi,-pi,-pi,-pi,-pi,-pi,-pi,-pi,-pi,-pi,-pi,-pi,-pi,-pi,0the lower bounds for the gridGRID_MAX=pi,pi,pi,pi,pi,pi,pi,pi,pi,pi,pi,pi,pi,pi,pi,pi,pi,pi,5the upper bounds for the gridGRID_WSTRIDE=5000frequency for dumping the gridWALKERS_MPI... PRINT( default=off ) Switch on MPI version of multiple walkers - not compatible with WALKERS_* options other than WALKERS_DIRFILE=COLVARthe name of the file on which to output these quantitiesARG=phi2,phi3,phi4,phi5,phi6,phi7,phi8,phi9,psi1,psi2,psi3,psi4,psi5,psi6,psi7,psi8,dihtrp_cacb,dihtrp_cbcg,ahc.lessthanthe input for this action is the scalar output from one or more other actions.STRIDE=200 PRINTcompulsory keyword ( default=1 ) the frequency with which the quantities of interest should be outputFILE=ENERGYthe name of the file on which to output these quantitiesARG=bias.bias,pb.biasthe input for this action is the scalar output from one or more other actions.STRIDE=200compulsory keyword ( default=1 ) the frequency with which the quantities of interest should be output
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 &
The former simulations should provide a converged (check for this) ensemble for the peptide. As shown in [73] 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. In this latter case we use a different noise model, i.e. NOISE=MOUTLIERS because the forward model for chemical shifts can result in systematic errors for some of them.
The following input file is meant to be appended to the former files.
# 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 MOLINFOPRINT ARG=jhanst.*,jhahnst.*,jw5ccyst.*,jw5ncyst.* STRIDE=2000 FILE=ST.JMOLTYPE=proteincompulsory keyword ( default=protein ) what kind of molecule is contained in the pdb file - usually not needed since protein/RNA/DNA are compatibleSTRUCTURE=egaawaass.pdb WHOLEMOLECULEScompulsory keyword a file in pdb format containing a reference structure.ENTITY0=1-111 # EEF1SB Implicit solvation #SETTINGS AUXFILE=user-doc/tutorials/others/isdb-1/reference-impl/index.ndx protein-h: GROUPthe atoms that make up a molecule that you wish to align.NDX_FILE=index.ndxthe name of index file (gromacs syntax)NDX_GROUP=Protein-H solv: EEFSOLVthe name of the group to be imported (gromacs syntax) - first group found is used by defaultATOMS=protein-hThe atoms to be included in the calculation, e.g.NOPBC( default=off ) ignore the periodic boundary conditions when calculating distancesNL_STRIDE=20compulsory keyword ( default=40 ) The frequency with which the neighbor list is updated.NL_BUFFER=0.1 bias: BIASVALUEcompulsory keyword ( default=0.1 ) The buffer to the intrinsic cutoff used when calculating pairwise interactions.ARG=solv # CVs, Psi9, Phi1 are not defined psi1: TORSIONthe input for this action is the scalar output from one or more other actions.ATOMS=@psi-1the four atoms involved in the torsional angleNOPBCpsi2: TORSION( default=off ) ignore the periodic boundary conditions when calculating distancesATOMS=@psi-2the four atoms involved in the torsional angleNOPBCpsi3: TORSION( default=off ) ignore the periodic boundary conditions when calculating distancesATOMS=@psi-3the four atoms involved in the torsional angleNOPBCpsi4: TORSION( default=off ) ignore the periodic boundary conditions when calculating distancesATOMS=@psi-4the four atoms involved in the torsional angleNOPBCpsi5: TORSION( default=off ) ignore the periodic boundary conditions when calculating distancesATOMS=@psi-5the four atoms involved in the torsional angleNOPBCpsi6: TORSION( default=off ) ignore the periodic boundary conditions when calculating distancesATOMS=@psi-6the four atoms involved in the torsional angleNOPBCpsi7: TORSION( default=off ) ignore the periodic boundary conditions when calculating distancesATOMS=@psi-7the four atoms involved in the torsional angleNOPBCpsi8: TORSION( default=off ) ignore the periodic boundary conditions when calculating distancesATOMS=@psi-8the four atoms involved in the torsional angleNOPBCphi2: TORSION( default=off ) ignore the periodic boundary conditions when calculating distancesATOMS=@phi-2the four atoms involved in the torsional angleNOPBCphi3: TORSION( default=off ) ignore the periodic boundary conditions when calculating distancesATOMS=@phi-3the four atoms involved in the torsional angleNOPBCphi4: TORSION( default=off ) ignore the periodic boundary conditions when calculating distancesATOMS=@phi-4the four atoms involved in the torsional angleNOPBCphi5: TORSION( default=off ) ignore the periodic boundary conditions when calculating distancesATOMS=@phi-5the four atoms involved in the torsional angleNOPBCphi6: TORSION( default=off ) ignore the periodic boundary conditions when calculating distancesATOMS=@phi-6the four atoms involved in the torsional angleNOPBCphi7: TORSION( default=off ) ignore the periodic boundary conditions when calculating distancesATOMS=@phi-7the four atoms involved in the torsional angleNOPBCphi8: TORSION( default=off ) ignore the periodic boundary conditions when calculating distancesATOMS=@phi-8the four atoms involved in the torsional angleNOPBCphi9: TORSION( default=off ) ignore the periodic boundary conditions when calculating distancesATOMS=@phi-9the four atoms involved in the torsional angleNOPBCahc: ALPHARMSD( default=off ) ignore the periodic boundary conditions when calculating distancesRESIDUES=allthis command is used to specify the set of residues that could conceivably form part of the secondary structure.TYPE=OPTIMALcompulsory keyword ( default=DRMSD ) the manner in which RMSD alignment is performed.LESS_THAN={RATIONAL R_0=0.12} # Bulky Trp residue dihedral dihtrp_cacb: TORSIONcalculate the number of variables less than a certain target value.ATOMS=67,47,49,52the four atoms involved in the torsional angleNOPBCdihtrp_cbcg: TORSION( default=off ) ignore the periodic boundary conditions when calculating distancesATOMS=47,49,52,53the four atoms involved in the torsional angleNOPBCprotein-ca: GROUP( default=off ) ignore the periodic boundary conditions when calculating distancesNDX_FILE=index.ndxthe name of index file (gromacs syntax)NDX_GROUP=C-alpha gyr: GYRATIONthe name of the group to be imported (gromacs syntax) - first group found is used by defaultTYPE=RADIUScompulsory keyword ( default=RADIUS ) The type of calculation relative to the Gyration Tensor you want to performATOMS=protein-cathe group of atoms that you are calculating the Gyration Tensor for.NOPBC# PBMetaD pb: PBMETAD ...( default=off ) ignore the periodic boundary conditions when calculating distancesARG=phi2,phi3,phi4,phi5,phi6,phi7,phi8,phi9,psi1,psi2,psi3,psi4,psi5,psi6,psi7,psi8,dihtrp_cacb,dihtrp_cbcg,ahc.lessthanthe input for this action is the scalar output from one or more other actions.SIGMA=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.01compulsory keyword the widths of the Gaussian hillsHEIGHT=0.5the height of the Gaussian hills, one for all biases.PACE=400compulsory keyword the frequency for hill addition, one for all biasesBIASFACTOR=20use well tempered metadynamics with this bias factor, one for all biases.GRID_MIN=-pi,-pi,-pi,-pi,-pi,-pi,-pi,-pi,-pi,-pi,-pi,-pi,-pi,-pi,-pi,-pi,-pi,-pi,0the lower bounds for the gridGRID_MAX=pi,pi,pi,pi,pi,pi,pi,pi,pi,pi,pi,pi,pi,pi,pi,pi,pi,pi,5the upper bounds for the gridGRID_WSTRIDE=5000frequency for dumping the gridWALKERS_MPI... PRINT( default=off ) Switch on MPI version of multiple walkers - not compatible with WALKERS_* options other than WALKERS_DIRFILE=COLVARthe name of the file on which to output these quantitiesARG=phi2,phi3,phi4,phi5,phi6,phi7,phi8,phi9,psi1,psi2,psi3,psi4,psi5,psi6,psi7,psi8,dihtrp_cacb,dihtrp_cbcg,ahc.lessthanthe input for this action is the scalar output from one or more other actions.STRIDE=200 PRINTcompulsory keyword ( default=1 ) the frequency with which the quantities of interest should be outputFILE=ENERGYthe name of the file on which to output these quantitiesARG=bias.bias,pb.biasthe input for this action is the scalar output from one or more other actions.STRIDE=200 # EXPERIMENTAL DATA SECTION # RDCs (Grzesiek et al.) # xGAAWAASS nh: RDC ...compulsory keyword ( default=1 ) the frequency with which the quantities of interest should be outputGYROM=-72.5388compulsory keyword ( default=1. ) Add the product of the gyromagnetic constants for the bond.SCALE=0.0001compulsory keyword ( default=1. ) Add the scaling factor to take into account concentration and other effects.NOPBC( default=off ) ignore the periodic boundary conditions when calculating distancesATOMS1=18,19the couple of atoms involved in each of the bonds for which you wish to calculate the RDC.COUPLING1=-5.4Add an experimental value for each coupling (needed by SVD and useful for \ref STATS).ATOMS2=25,26the couple of atoms involved in each of the bonds for which you wish to calculate the RDC.COUPLING2=-1.26Add an experimental value for each coupling (needed by SVD and useful for \ref STATS).ATOMS3=35,36the couple of atoms involved in each of the bonds for which you wish to calculate the RDC.COUPLING3=-5.22Add an experimental value for each coupling (needed by SVD and useful for \ref STATS).ATOMS4=45,46the couple of atoms involved in each of the bonds for which you wish to calculate the RDC.COUPLING4=-0.91Add an experimental value for each coupling (needed by SVD and useful for \ref STATS).ATOMS5=69,70the couple of atoms involved in each of the bonds for which you wish to calculate the RDC.COUPLING5=2.33Add an experimental value for each coupling (needed by SVD and useful for \ref STATS).ATOMS6=79,80the couple of atoms involved in each of the bonds for which you wish to calculate the RDC.COUPLING6=-2.88Add an experimental value for each coupling (needed by SVD and useful for \ref STATS).ATOMS7=89,90the couple of atoms involved in each of the bonds for which you wish to calculate the RDC.COUPLING7=-8.37Add an experimental value for each coupling (needed by SVD and useful for \ref STATS).ATOMS8=100,101the couple of atoms involved in each of the bonds for which you wish to calculate the RDC.COUPLING8=-3.78 ... # ExAAWAASx caha: RDC ...Add an experimental value for each coupling (needed by SVD and useful for \ref STATS).GYROM=179.9319compulsory keyword ( default=1. ) Add the product of the gyromagnetic constants for the bond.SCALE=0.0001compulsory keyword ( default=1. ) Add the scaling factor to take into account concentration and other effects.NOPBC( default=off ) ignore the periodic boundary conditions when calculating distancesATOMS1=5,6the couple of atoms involved in each of the bonds for which you wish to calculate the RDC.COUPLING1=12.95Add an experimental value for each coupling (needed by SVD and useful for \ref STATS).ATOMS2=27,28the couple of atoms involved in each of the bonds for which you wish to calculate the RDC.COUPLING2=11.5Add an experimental value for each coupling (needed by SVD and useful for \ref STATS).ATOMS3=37,38the couple of atoms involved in each of the bonds for which you wish to calculate the RDC.COUPLING3=21.42Add an experimental value for each coupling (needed by SVD and useful for \ref STATS).ATOMS4=47,48the couple of atoms involved in each of the bonds for which you wish to calculate the RDC.COUPLING4=-9.37Add an experimental value for each coupling (needed by SVD and useful for \ref STATS).ATOMS5=71,72the couple of atoms involved in each of the bonds for which you wish to calculate the RDC.COUPLING5=10.01Add an experimental value for each coupling (needed by SVD and useful for \ref STATS).ATOMS6=81,82the couple of atoms involved in each of the bonds for which you wish to calculate the RDC.COUPLING6=15.01Add an experimental value for each coupling (needed by SVD and useful for \ref STATS).ATOMS7=91,92the couple of atoms involved in each of the bonds for which you wish to calculate the RDC.COUPLING7=15.73 ... #RDCS byrdcnh: METAINFERENCE ...Add an experimental value for each coupling (needed by SVD and useful for \ref STATS).ARG=(nh\.rdc-.*),pb.biasthe input for this action is the scalar output from one or more other actions.PARARG=(nh\.exp-.*)reference values for the experimental data, these can be provided as arguments without derivativesREWEIGHT( default=off ) simple REWEIGHT using the latest ARG as energyNOISETYPE=MGAUSScompulsory keyword ( default=MGAUSS ) functional form of the noise (GAUSS,MGAUSS,OUTLIERS,MOUTLIERS,GENERIC)OPTSIGMAMEAN=SEM_MAXcompulsory keyword ( default=NONE ) Set to NONE/SEM to manually set sigma mean, or to estimate it on the flyAVERAGING=400Stride for calculation of averaged weights and sigma_meanSCALEDATA( default=off ) Set to TRUE if you want to sample a scaling factor common to all values and replicasSCALE_PRIOR=GAUSSIANcompulsory keyword ( default=FLAT ) either FLAT or GAUSSIANSCALE0=8.0could not find this keywordDSCALE=0.5maximum MC move of the scaling factorSIGMA0=5.0could not find this keywordSIGMA_MIN=0.0001compulsory keyword ( default=0.0 ) minimum value of the uncertainty parameterSIGMA_MAX=15.0compulsory keyword ( default=10. ) maximum value of the uncertainty parameterDSIGMA=0.1maximum MC move of the uncertainty parameterWRITE_STRIDE=10000 ... #RDCS byrdccaha: METAINFERENCE ...compulsory keyword ( default=10000 ) write the status to a file every N steps, this can be used for restart/continuationARG=(caha\.rdc-.*),pb.biasthe input for this action is the scalar output from one or more other actions.PARARG=(caha\.exp-.*)reference values for the experimental data, these can be provided as arguments without derivativesREWEIGHT( default=off ) simple REWEIGHT using the latest ARG as energyNOISETYPE=MGAUSScompulsory keyword ( default=MGAUSS ) functional form of the noise (GAUSS,MGAUSS,OUTLIERS,MOUTLIERS,GENERIC)OPTSIGMAMEAN=SEM_MAXcompulsory keyword ( default=NONE ) Set to NONE/SEM to manually set sigma mean, or to estimate it on the flyAVERAGING=400Stride for calculation of averaged weights and sigma_meanSCALEDATA( default=off ) Set to TRUE if you want to sample a scaling factor common to all values and replicasSCALE_PRIOR=GAUSSIANcompulsory keyword ( default=FLAT ) either FLAT or GAUSSIANSCALE0=9.0could not find this keywordDSCALE=0.5maximum MC move of the scaling factorSIGMA0=5.0could not find this keywordSIGMA_MIN=0.0001compulsory keyword ( default=0.0 ) minimum value of the uncertainty parameterSIGMA_MAX=15.0compulsory keyword ( default=10. ) maximum value of the uncertainty parameterDSIGMA=0.1maximum MC move of the uncertainty parameterWRITE_STRIDE=10000 ... # xGxAWxASx jhan: JCOUPLING ...compulsory keyword ( default=10000 ) write the status to a file every N steps, this can be used for restart/continuationTYPE=HANcompulsory keyword Type of J-coupling to compute (HAN,HAHN,CCG,NCG,CUSTOM)NOPBC( default=off ) ignore the periodic boundary conditions when calculating distancesATOMS1=@psi-2the 4 atoms involved in each of the bonds for which you wish to calculate the J-coupling.COUPLING1=-0.49Add an experimental value for each coupling You can use multiple instances of this keyword i.e.ATOMS2=@psi-4the 4 atoms involved in each of the bonds for which you wish to calculate the J-coupling.COUPLING2=-0.54Add an experimental value for each coupling You can use multiple instances of this keyword i.e.ATOMS3=@psi-5the 4 atoms involved in each of the bonds for which you wish to calculate the J-coupling.COUPLING3=-0.53Add an experimental value for each coupling You can use multiple instances of this keyword i.e.ATOMS4=@psi-7the 4 atoms involved in each of the bonds for which you wish to calculate the J-coupling.COUPLING4=-0.39Add an experimental value for each coupling You can use multiple instances of this keyword i.e.ATOMS5=@psi-8the 4 atoms involved in each of the bonds for which you wish to calculate the J-coupling.COUPLING5=-0.39 ... # xxAAWAASS jhahn: JCOUPLING ...Add an experimental value for each coupling You can use multiple instances of this keyword i.e.TYPE=HAHNcompulsory keyword Type of J-coupling to compute (HAN,HAHN,CCG,NCG,CUSTOM)NOPBC( default=off ) ignore the periodic boundary conditions when calculating distancesATOMS1=@phi-2the 4 atoms involved in each of the bonds for which you wish to calculate the J-coupling.COUPLING1=6.05Add an experimental value for each coupling You can use multiple instances of this keyword i.e.ATOMS2=@phi-3the 4 atoms involved in each of the bonds for which you wish to calculate the J-coupling.COUPLING2=5.95Add an experimental value for each coupling You can use multiple instances of this keyword i.e.ATOMS3=@phi-4the 4 atoms involved in each of the bonds for which you wish to calculate the J-coupling.COUPLING3=6.44Add an experimental value for each coupling You can use multiple instances of this keyword i.e.ATOMS4=@phi-5the 4 atoms involved in each of the bonds for which you wish to calculate the J-coupling.COUPLING4=6.53Add an experimental value for each coupling You can use multiple instances of this keyword i.e.ATOMS5=@phi-6the 4 atoms involved in each of the bonds for which you wish to calculate the J-coupling.COUPLING5=5.93Add an experimental value for each coupling You can use multiple instances of this keyword i.e.ATOMS6=@phi-7the 4 atoms involved in each of the bonds for which you wish to calculate the J-coupling.COUPLING6=6.98Add an experimental value for each coupling You can use multiple instances of this keyword i.e.ATOMS7=@phi-8the 4 atoms involved in each of the bonds for which you wish to calculate the J-coupling.COUPLING7=7.16 ... # xxxxWxxxx jccg: JCOUPLING ...Add an experimental value for each coupling You can use multiple instances of this keyword i.e.TYPE=CCGcompulsory keyword Type of J-coupling to compute (HAN,HAHN,CCG,NCG,CUSTOM)NOPBC( default=off ) ignore the periodic boundary conditions when calculating distancesATOMS1=@chi1-5the 4 atoms involved in each of the bonds for which you wish to calculate the J-coupling.COUPLING1=1.59 ... # xxxxWxxxx jncg: JCOUPLING ...Add an experimental value for each coupling You can use multiple instances of this keyword i.e.TYPE=NCGcompulsory keyword Type of J-coupling to compute (HAN,HAHN,CCG,NCG,CUSTOM)NOPBC( default=off ) ignore the periodic boundary conditions when calculating distancesATOMS1=@chi1-5the 4 atoms involved in each of the bonds for which you wish to calculate the J-coupling.COUPLING1=1.21 ... #JC METAINFERENCE byj: METAINFERENCE ...Add an experimental value for each coupling You can use multiple instances of this keyword i.e.ARG=(jhan\.j-.*),(jhahn\.j-.*),(jccg\.j.*),(jncg\.j.*),pb.biasthe input for this action is the scalar output from one or more other actions.PARARG=(jhan\.exp-.*),(jhahn\.exp-.*),(jccg\.exp.*),(jncg\.exp.*)reference values for the experimental data, these can be provided as arguments without derivativesREWEIGHT( default=off ) simple REWEIGHT using the latest ARG as energyNOISETYPE=MGAUSScompulsory keyword ( default=MGAUSS ) functional form of the noise (GAUSS,MGAUSS,OUTLIERS,MOUTLIERS,GENERIC)OPTSIGMAMEAN=SEM_MAXcompulsory keyword ( default=NONE ) Set to NONE/SEM to manually set sigma mean, or to estimate it on the flyAVERAGING=400Stride for calculation of averaged weights and sigma_meanSIGMA0=5.0could not find this keywordSIGMA_MIN=0.0001compulsory keyword ( default=0.0 ) minimum value of the uncertainty parameterSIGMA_MAX=15.0compulsory keyword ( default=10. ) maximum value of the uncertainty parameterDSIGMA=0.1maximum MC move of the uncertainty parameterWRITE_STRIDE=10000 ... # Chemical shifts #SETTINGS AUXFOLDER=user-doc/tutorials/others/isdb-1/m_and_m/data # Chemical shifts cs: CS2BACKBONE ...compulsory keyword ( default=10000 ) write the status to a file every N steps, this can be used for restart/continuationNOPBC( default=off ) ignore the periodic boundary conditions when calculating distancesATOMS=1-111The atoms to be included in the calculation, e.g.DATADIR=datacompulsory keyword ( default=data/ ) The folder with the experimental chemical shifts.TEMPLATE=egaawaass.pdbcompulsory keyword ( default=template.pdb ) A PDB file of the protein system.DOSCORE( default=off ) activate metainferenceARG=pb.biasthe input for this action is the scalar output from one or more other actions.NOISETYPE=MGAUSScompulsory keyword ( default=MGAUSS ) functional form of the noise (GAUSS,MGAUSS,OUTLIERS,MOUTLIERS,GENERIC)REWEIGHT( default=off ) simple REWEIGHT using the ARG as energyOPTSIGMAMEAN=SEM_MAXcompulsory keyword ( default=NONE ) Set to NONE/SEM to manually set sigma mean, or to estimate it on the flyAVERAGING=400Stride for calculation of averaged weights and sigma_meanSIGMA0=4.0could not find this keywordSIGMA_MIN=0.0001compulsory keyword ( default=0.0 ) minimum value of the uncertainty parameterSIGMA_MAX=5.00compulsory keyword ( default=10. ) maximum value of the uncertainty parameterWRITE_STRIDE=10000 ... mcs: BIASVALUEcompulsory keyword ( default=10000 ) write the status to a file every N steps, this can be used for restart/continuationARG=cs.score # output from METAINFERENCE PRINTthe input for this action is the scalar output from one or more other actions.ARG=byrdcnh.*the input for this action is the scalar output from one or more other actions.STRIDE=200compulsory keyword ( default=1 ) the frequency with which the quantities of interest should be outputFILE=BAYES.RDC.NH PRINTthe name of the file on which to output these quantitiesARG=byrdccaha.*the input for this action is the scalar output from one or more other actions.STRIDE=200compulsory keyword ( default=1 ) the frequency with which the quantities of interest should be outputFILE=BAYES.RDC.CAHA PRINTthe name of the file on which to output these quantitiesARG=byj.*the input for this action is the scalar output from one or more other actions.STRIDE=200compulsory keyword ( default=1 ) the frequency with which the quantities of interest should be outputFILE=BAYES.J PRINTthe name of the file on which to output these quantitiesARG=cs.*the input for this action is the scalar output from one or more other actions.STRIDE=200compulsory keyword ( default=1 ) the frequency with which the quantities of interest should be outputFILE=BAYES.CS # the following are useful for the analysis on-the-fly of the quality of the agreement with the experimentl data ens: ENSEMBLE ...the name of the file on which to output these quantitiesARG=(nh\.rdc-.*),(caha\.rdc-.*),(jhan\.j-.*),(jhahn\.j-.*),(jccg\.j-.*),(jncg\.j-.*),(cs\...-.*),pb.biasthe input for this action is the scalar output from one or more other actions.REWEIGHT... nhst: STATS ...( default=off ) simple REWEIGHT using the latest ARG as energyARG=(ens\.nh\.rdc-.*)the input for this action is the scalar output from one or more other actions.PARARG=(nh\.exp-.*) ... cahast: STATS ...the input for this action is the scalar output from one or more other actions without derivatives.ARG=(ens\.caha\.rdc-.*)the input for this action is the scalar output from one or more other actions.PARARG=(caha\.exp-.*) ... csst: STATS ...the input for this action is the scalar output from one or more other actions without derivatives.ARG=(ens\.cs\...-.*)the input for this action is the scalar output from one or more other actions.PARARG=(cs\.exp.*) ... jhanst: STATS ...the input for this action is the scalar output from one or more other actions without derivatives.ARG=(ens\.jhan\.j-.*)the input for this action is the scalar output from one or more other actions.PARARG=(jhan\.exp-.*) ... jhahnst: STATS ...the input for this action is the scalar output from one or more other actions without derivatives.ARG=(ens\.jhahn\.j-.*)the input for this action is the scalar output from one or more other actions.PARARG=(jhahn\.exp-.*) ... jw5ccyst: STATS ...the input for this action is the scalar output from one or more other actions without derivatives.ARG=(ens\.jccg\.j.*),(ens\.jccg\.j.*)the input for this action is the scalar output from one or more other actions.PARARG=(jccg\.exp-.*),(jccg\.exp-.*)the input for this action is the scalar output from one or more other actions without derivatives.SQDEVSUM... jw5ncyst: STATS ...( default=off ) calculates only SQDEVSUMARG=(ens\.jncg\.j.*),(ens\.jncg\.j.*)the input for this action is the scalar output from one or more other actions.PARARG=(jncg\.exp-.*),(jncg\.exp-.*)the input for this action is the scalar output from one or more other actions without derivatives.SQDEVSUM... #output from STATS PRINT( default=off ) calculates only SQDEVSUMARG=nhst.*the input for this action is the scalar output from one or more other actions.STRIDE=2000compulsory keyword ( default=1 ) the frequency with which the quantities of interest should be outputFILE=ST.RDC.NH PRINTthe name of the file on which to output these quantitiesARG=cahast.*the input for this action is the scalar output from one or more other actions.STRIDE=2000compulsory keyword ( default=1 ) the frequency with which the quantities of interest should be outputFILE=ST.RDC.CAHA PRINTthe name of the file on which to output these quantitiesARG=csst.*the input for this action is the scalar output from one or more other actions.STRIDE=2000compulsory keyword ( default=1 ) the frequency with which the quantities of interest should be outputFILE=ST.CSthe name of the file on which to output these quantities
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 [19] .
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 &