Inconsistency in protein-peptide modeling with CHARMM

GROMACS version: 2021.3
GROMACS modification: Yes/No

I am trying to do modeling of a small peptide with a protein. I tried two ways of treating the complex:

  1. small peptide as a part of a protein
  2. small peptide as a ligand

In the first case I don’t have to do preprocessing and I can just start the calculation as for one simple protein. In second case I am obtaining the .itp files for a small peptide from SwissParam and start calculation as a protein-ligand complex.

I expected to see both results to be similar, however every time I treat peptide as a ligand (and not part of a protein) , the peptide immediately looses its secondary structure, got unfolded and produce results far away from crystal structure.

Why there is such a discrepancy and what way of treating peptide is more correct?


The peptide topology should be generated using the standard protein force field, not SwissParam. The latter is designed for general-purpose ligand topologies and is not highly optimized in the way that the protein force field requires. pdb2gmx can generate the entire topology for you.

@jalemkul, thank you for the explanation! My problem however that I have non-standard peptides with the bridges between some of the residues, so I can not treat them as a part of a protein. Maybe it will help if I increase the penalty in POSRES_LIGAND file?

A low quality topology will eventually cause problems. You should parametrize the only nonstandard residues carefully and incorporate them into the protein force field so you can generate a reasonable topology.

@jalemkul , thanks. I have red a brief explanation on how to create a forcefield for a protein with some non-standard residues, but my peptides have a ring system (please, see attached), so what could be the workflow in this case? bridgepept

Parametrize the individual nonstandard amino acids, and add an entry in specbond.dat that specifies how they are bonded, and pdb2gmx will create the connection between the two.