Martini 2.2 topologies for alpha-helical peptides _with uncharged termini_ are incorrect

[Ladme]

Context This issue is related to issues #560 and #566 (yes, me again). Thanks to your efforts, it seems that Martini 2.2 topologies of peptides with charged termini are now generated correctly. However, alpha-helical peptides with uncharged termini are assigned incorrect bead types for their terminal backbone beads.

Description The current dev version of vermouth-martinize (version 0.10.1.dev6) does not assign correct backbone beads to terminal amino acids when these amino acids are designated as being uncharged using the -nt flag.

Using martinize2 -f l20.pdb -o topol.top -x l20_martini.pdb -ss HHHHHHHHHHHHHHHHHHHH -p backbone -ff martini22 -nt (where l20.pdb is a MODELLER-constructed 20 leucines long alpha-helical peptide), generates the following [ atoms ] block:

[ atoms ]
 1 P5   1 LEU BB   1   0
 2 AC1  1 LEU SC1  2 0.0
 3 Nd   2 LEU BB   3 0.0
 4 AC1  2 LEU SC1  4 0.0
 5 Nd   3 LEU BB   5 0.0
 6 AC1  3 LEU SC1  6 0.0
 7 Nd   4 LEU BB   7 0.0
 8 AC1  4 LEU SC1  8 0.0
 9 N0   5 LEU BB   9 0.0
10 AC1  5 LEU SC1 10 0.0
11 N0   6 LEU BB  11 0.0
12 AC1  6 LEU SC1 12 0.0
13 N0   7 LEU BB  13 0.0
14 AC1  7 LEU SC1 14 0.0
15 N0   8 LEU BB  15 0.0
16 AC1  8 LEU SC1 16 0.0
17 N0   9 LEU BB  17 0.0
18 AC1  9 LEU SC1 18 0.0
19 N0  10 LEU BB  19 0.0
20 AC1 10 LEU SC1 20 0.0
21 N0  11 LEU BB  21 0.0
22 AC1 11 LEU SC1 22 0.0
23 N0  12 LEU BB  23 0.0
24 AC1 12 LEU SC1 24 0.0
25 N0  13 LEU BB  25 0.0
26 AC1 13 LEU SC1 26 0.0
27 N0  14 LEU BB  27 0.0
28 AC1 14 LEU SC1 28 0.0
29 N0  15 LEU BB  29 0.0
30 AC1 15 LEU SC1 30 0.0
31 N0  16 LEU BB  31 0.0
32 AC1 16 LEU SC1 32 0.0
33 Na  17 LEU BB  33 0.0
34 AC1 17 LEU SC1 34 0.0
35 Na  18 LEU BB  35 0.0
36 AC1 18 LEU SC1 36 0.0
37 Na  19 LEU BB  37 0.0
38 AC1 19 LEU SC1 38 0.0
39 P5  20 LEU BB  39   0
40 AC1 20 LEU SC1 40 0.0

Meanwhile, using martinize.py script version 2.6_3 (python martinize.py -f l20.pdb -o topol.top -x l20_martini.pdb -ss HHHHHHHHHHHHHHHHHHHH -p backbone -nt) generates this [ atom ] block:

[ atoms ]
    1    Nd     1   LEU    BB     1  0.0000 ; 1
    2   AC1     1   LEU   SC1     2  0.0000 ; 1
    3    Nd     2   LEU    BB     3  0.0000 ; 1
    4   AC1     2   LEU   SC1     4  0.0000 ; 1
    5    Nd     3   LEU    BB     5  0.0000 ; 1
    6   AC1     3   LEU   SC1     6  0.0000 ; 1
    7    Nd     4   LEU    BB     7  0.0000 ; 1
    8   AC1     4   LEU   SC1     8  0.0000 ; 1
    9    N0     5   LEU    BB     9  0.0000 ; H
   10   AC1     5   LEU   SC1    10  0.0000 ; H
   11    N0     6   LEU    BB    11  0.0000 ; H
   12   AC1     6   LEU   SC1    12  0.0000 ; H
   13    N0     7   LEU    BB    13  0.0000 ; H
   14   AC1     7   LEU   SC1    14  0.0000 ; H
   15    N0     8   LEU    BB    15  0.0000 ; H
   16   AC1     8   LEU   SC1    16  0.0000 ; H
   17    N0     9   LEU    BB    17  0.0000 ; H
   18   AC1     9   LEU   SC1    18  0.0000 ; H
   19    N0    10   LEU    BB    19  0.0000 ; H
   20   AC1    10   LEU   SC1    20  0.0000 ; H
   21    N0    11   LEU    BB    21  0.0000 ; H
   22   AC1    11   LEU   SC1    22  0.0000 ; H
   23    N0    12   LEU    BB    23  0.0000 ; H
   24   AC1    12   LEU   SC1    24  0.0000 ; H
   25    N0    13   LEU    BB    25  0.0000 ; H
   26   AC1    13   LEU   SC1    26  0.0000 ; H
   27    N0    14   LEU    BB    27  0.0000 ; H
   28   AC1    14   LEU   SC1    28  0.0000 ; H
   29    N0    15   LEU    BB    29  0.0000 ; H
   30   AC1    15   LEU   SC1    30  0.0000 ; H
   31    N0    16   LEU    BB    31  0.0000 ; H
   32   AC1    16   LEU   SC1    32  0.0000 ; H
   33    Na    17   LEU    BB    33  0.0000 ; 2
   34   AC1    17   LEU   SC1    34  0.0000 ; 2
   35    Na    18   LEU    BB    35  0.0000 ; 2
   36   AC1    18   LEU   SC1    36  0.0000 ; 2
   37    Na    19   LEU    BB    37  0.0000 ; 2
   38   AC1    19   LEU   SC1    38  0.0000 ; 2
   39    Na    20   LEU    BB    39  0.0000 ; 2
   40   AC1    20   LEU   SC1    40  0.0000 ; 2

There is a discrepancy between the type of beads number 1 and 39 which should be Nd and Na, respectively, not P5 as in vermouth-martinize.

[pckroon]

Hello hello! Always happy with (detailed) bug reports! It means that 1) people are using martinize2, and 2) actually check their topologies, and 3) care enough about martinize2 to ask for it to be fixed :)

Alright, this is indeed related to #560 (more than #566). However, this issue is the other way around. In #560 the issue was that Links (the sec-struct specific parameters) would override the terminal modifications, while here the issue is that the neutral terminal modification overrides the Link parameters. So I'm going to ask the same question as before: what's the correct behaviour? I have the feeling that the P5 is actually the correct answer, the terminal backbone should be way more polar than internal backbones, no? So, paging @fgrunewald @paulocts

[pckroon]

Hello hello! I had some mail contact with Paulo, and the verdict is in: you found a bug in martinize 1! I'll copy-paste part of his mail here:

Unfortunately, the original papers are not very clear. But we still can see the original code (before martinize1). Take a look here if you are brave enough: http://cgmartini.nl/images/tools/seq2itp/seq2itp_martini2.1.pl. I may be wrong, but what is defined in this code as BEAD_BBN_BND is the neutral termini ... and it is P5. But I am not sure. If we think about what would be chemically expected, I also think it should be P5. In addition, I think that the behavior of TM peptides is better described with P5 as termini as it helps stabilize them as TM. So, there is a pragmatic reason to also keep as P5.

Also note that a protein with a helix terminus is not very physically relevant: DSSP tends to assign the termini C, and the ends of helices will partially unfold.