Подсчет количества водородных связей нескольких файлов PDB

Я пытался найти систематический способ подсчета количества водородных связей для нескольких файлов PDB.

DSSP показывает мне общее количество водородных связей, когда я делаю dssp.exe 1A11.pdb, давая мне общее количество 24(я прав, хорошо?)

24 96.0   TOTAL NUMBER OF HYDROGEN BONDS OF TYPE O(I)-->H-N(J)  , SAME NUMBER PER 100 RESIDUES

Но когда я запускаю эту команду в R:

> require(biod3d)
> x=dssp(read.pdb("1A11"), full=T)

у него не похоже количество водородных связей. Это показывает следующее:

> x$hbonds
     BP1 BP2 NH-O.1   E1 O-HN.1   E2 NH-O.2   E3 O-HN.2   E4 ChainBP1 ChainBP2 Chain1 Chain2 Chain3 Chain4
1-A   NA  NA     NA  0.0      5 -2.1     NA  0.0      4 -0.3     <NA>     <NA>   <NA>      A   <NA>      A
2-A   NA  NA      3 -0.2      6 -1.5      4 -0.2      7 -0.2     <NA>     <NA>      A      A      A      A
3-A   NA  NA      4 -0.2      7 -1.8      5 -0.2      2 -0.2     <NA>     <NA>      A      A      A      A
4-A   NA  NA      1 -0.3      8 -1.5      6 -0.2      7 -0.6     <NA>     <NA>      A      A      A      A
5-A   NA  NA      1 -2.1      9 -1.9      6 -0.2     10 -0.2     <NA>     <NA>      A      A      A      A
6-A   NA  NA      2 -1.5     10 -1.3      1 -0.3      5 -0.2     <NA>     <NA>      A      A      A      A
7-A   NA  NA      3 -1.8     11 -1.3      4 -0.6      5 -0.2     <NA>     <NA>      A      A      A      A
8-A   NA  NA      4 -1.5     12 -1.2      3 -0.2      6 -0.2     <NA>     <NA>      A      A      A      A
9-A   NA  NA      5 -1.9     13 -1.1      4 -0.2      7 -0.2     <NA>     <NA>      A      A      A      A
10-A  NA  NA      6 -1.3     14 -2.3      5 -0.2     13 -0.4     <NA>     <NA>      A      A      A      A
11-A  NA  NA      7 -1.3     15 -1.8     12 -0.2     16 -0.3     <NA>     <NA>      A      A      A      A
12-A  NA  NA      8 -1.2     16 -1.8     13 -0.2     11 -0.2     <NA>     <NA>      A      A      A      A
13-A  NA  NA      9 -1.1     17 -2.4     10 -0.4     18 -0.5     <NA>     <NA>      A      A      A      A
14-A  NA  NA     10 -2.3     18 -1.6     15 -0.2     13 -0.2     <NA>     <NA>      A      A      A      A
15-A  NA  NA     11 -1.8     19 -1.8     17 -0.2     20 -0.5     <NA>     <NA>      A      A      A      A
16-A  NA  NA     12 -1.8     20 -0.8     17 -0.3     14 -0.2     <NA>     <NA>      A      A      A      A
17-A  NA  NA     13 -2.4     21 -1.9     19 -0.2     16 -0.3     <NA>     <NA>      A      A      A      A
18-A  NA  NA     14 -1.6     22 -2.0     13 -0.5     21 -1.0     <NA>     <NA>      A      A      A      A
19-A  NA  NA     15 -1.8     23 -2.0     20 -0.3     24 -0.3     <NA>     <NA>      A      A      A      A
20-A  NA  NA     16 -0.8     24 -1.5     15 -0.5     19 -0.3     <NA>     <NA>      A      A      A      A
21-A  NA  NA     17 -1.9     25 -1.0     18 -1.0     19 -0.2     <NA>     <NA>      A      A      A      A
22-A  NA  NA     18 -2.0     25 -0.3     23 -0.2     20 -0.2     <NA>     <NA>      A      A      A      A
23-A  NA  NA     19 -2.0     22 -0.2     18 -0.3     21 -0.2     <NA>     <NA>      A      A      A      A
24-A  NA  NA     20 -1.5     23 -0.2     19 -0.3     22 -0.2     <NA>     <NA>      A      A      A      A
25-A  NA  NA     21 -1.0     23 -0.2     22 -0.3     24 -0.2     <NA>     <NA>      A      A      A      A

Итак, мой вопрос: общее количество водородных связей скрыто где-то в этом выводе или где-то еще (или не существует)?

Если мне удастся подсчитать количество водородных связей в одной PDB с помощью Bio3D и R, я смогу создать цикл и перебрать все оставшиеся файлы PDB.

Спасибо за все!

Ответы (2)

Число водородных связей на самом деле не может быть определено только с помощью программного обеспечения. Но вы можете попробовать использовать STRIDE ( http://structure.usc.edu/stride/ ), который берет имя файла с флагом -h для вывода количества водородных связей. Затем вы можете написать небольшой сценарий оболочки, который будет передавать каждый файл и хранить нужные данные в любом удобном для вас формате.

Это правда! Я думал, что может существовать прямой способ сделать это :)

Это первоначальное решение, но его можно улучшить:

out = system2("f:/stride/stride.exe", "f:/stride/1a11.pdb -h", stdout=T)
nh  = strsplit(out[grepl("HBT", out)], "  ")[[1]][2]
print(nh) # Output: 26 (hydrogen bonds)

Я скомпилировал версию Stride для Windows и загрузил сюда: http://lcrserver.icmc.usp.br/~daniel/bin/stride.exe