Как изменить содержимое столбцов в файле PDB белка по ключам и значениям, хранящимся в текстовом файле, используя только инструменты bash

Question

как написано в заголовке этого поста, я пытаюсь заменить значения в так называемых столбцах B (температурный коэффициент) и / или q (заполненность) файла PDB (банк данных белка), в котором хранятся трехмерные координатыкаждый атом, принадлежащий этому белку.Простой CSV-файл с двумя столбцами - это исходный файл, содержащий такие значения (во втором столбце) для конкретной аминокислоты (порядковый номер в первом столбце)

Краткий пример ограниченного исходного файла (source.csv)к первым двум аминокислотам (в реальном наборе данных их сотни):

1, 11.25
2, 16.49

Соответствующий блок файла назначения (dest.pdb):

ATOM      1  N   MET A   1     105.382 119.360 102.631  1.00   0.00
ATOM      2 CA   MET A   1     105.155 118.751 103.942  1.00   0.00
ATOM      3 HA   MET A   1     104.645 119.496 104.551  1.00   0.00
ATOM      4 CB   MET A   1     104.212 117.542 103.804  1.00   0.00
ATOM      5HB1   MET A   1     104.120 117.057 104.775  1.00   0.00
ATOM      6HB2   MET A   1     104.631 116.826 103.095  1.00   0.00
ATOM      7 CG   MET A   1     102.801 117.937 103.353  1.00   0.00
ATOM      8HG1   MET A   1     102.862 118.327 102.336  1.00   0.00
ATOM      9HG2   MET A   1     102.436 118.736 103.999  1.00   0.00
ATOM     10 SD   MET A   1     101.579 116.590 103.371  1.00   0.00
ATOM     11 CE   MET A   1     101.404 116.275 105.156  1.00   0.00
ATOM     12HE1   MET A   1     100.603 115.555 105.325  1.00   0.00
ATOM     13HE2   MET A   1     102.327 115.865 105.565  1.00   0.00
ATOM     14HE3   MET A   1     101.158 117.201 105.676  1.00   0.00
ATOM     15  C   MET A   1     106.423 118.387 104.697  1.00   0.00
ATOM     16  O   MET A   1     107.511 118.334 104.134  1.00   0.00
ATOM     17  N   GLU A   2     106.296 118.095 105.999  1.00   0.00
ATOM     18  H   GLU A   2     105.398 118.148 106.454  1.00   0.00
ATOM     19 CA   GLU A   2     107.495 117.802 106.786  1.00   0.00
ATOM     20 HA   GLU A   2     108.068 118.718 106.664  1.00   0.00
ATOM     21 CB   GLU A   2     107.242 117.714 108.295  1.00   0.00
ATOM     22HB1   GLU A   2     106.839 116.732 108.520  1.00   0.00
ATOM     23HB2   GLU A   2     106.494 118.455 108.581  1.00   0.00
ATOM     24 CG   GLU A   2     108.519 117.970 109.128  1.00   0.00
ATOM     25HG1   GLU A   2     108.323 117.660 110.155  1.00   0.00
ATOM     26HG2   GLU A   2     109.328 117.336 108.762  1.00   0.00
ATOM     27 CD   GLU A   2     109.002 119.432 109.126  1.00   0.00
ATOM     28OE1   GLU A   2     109.449 119.916 108.058  1.00   0.00
ATOM     29OE2   GLU A   2     109.026 120.057 110.206  1.00   0.00
ATOM     30  C   GLU A   2     108.446 116.757 106.163  1.00   0.00
ATOM     31  O   GLU A   2     109.650 117.015 106.154  1.00   0.00

Мне нужночтобы получить файл результатов (result.pdb):

ATOM      1  N   MET A   1     105.382 119.360 102.631  1.00  11.25
ATOM      2 CA   MET A   1     105.155 118.751 103.942  1.00  11.25
ATOM      3 HA   MET A   1     104.645 119.496 104.551  1.00  11.25
ATOM      4 CB   MET A   1     104.212 117.542 103.804  1.00  11.25
ATOM      5HB1   MET A   1     104.120 117.057 104.775  1.00  11.25
ATOM      6HB2   MET A   1     104.631 116.826 103.095  1.00  11.25
ATOM      7 CG   MET A   1     102.801 117.937 103.353  1.00  11.25
ATOM      8HG1   MET A   1     102.862 118.327 102.336  1.00  11.25
ATOM      9HG2   MET A   1     102.436 118.736 103.999  1.00  11.25
ATOM     10 SD   MET A   1     101.579 116.590 103.371  1.00  11.25
ATOM     11 CE   MET A   1     101.404 116.275 105.156  1.00  11.25
ATOM     12HE1   MET A   1     100.603 115.555 105.325  1.00  11.25
ATOM     13HE2   MET A   1     102.327 115.865 105.565  1.00  11.25
ATOM     14HE3   MET A   1     101.158 117.201 105.676  1.00  11.25
ATOM     15  C   MET A   1     106.423 118.387 104.697  1.00  11.25
ATOM     16  O   MET A   1     107.511 118.334 104.134  1.00  11.25
ATOM     17  N   GLU A   2     106.296 118.095 105.999  1.00  16.49
ATOM     18  H   GLU A   2     105.398 118.148 106.454  1.00  16.49
ATOM     19 CA   GLU A   2     107.495 117.802 106.786  1.00  16.49
ATOM     20 HA   GLU A   2     108.068 118.718 106.664  1.00  16.49
ATOM     21 CB   GLU A   2     107.242 117.714 108.295  1.00  16.49
ATOM     22HB1   GLU A   2     106.839 116.732 108.520  1.00  16.49
ATOM     23HB2   GLU A   2     106.494 118.455 108.581  1.00  16.49
ATOM     24 CG   GLU A   2     108.519 117.970 109.128  1.00  16.49
ATOM     25HG1   GLU A   2     108.323 117.660 110.155  1.00  16.49
ATOM     26HG2   GLU A   2     109.328 117.336 108.762  1.00  16.49
ATOM     27 CD   GLU A   2     109.002 119.432 109.126  1.00  16.49
ATOM     28OE1   GLU A   2     109.449 119.916 108.058  1.00  16.49
ATOM     29OE2   GLU A   2     109.026 120.057 110.206  1.00  16.49
ATOM     30  C   GLU A   2     108.446 116.757 106.163  1.00  16.49
ATOM     31  O   GLU A   2     109.650 117.015 106.154  1.00  16.49

, в котором оценки аминокислоты применяются ко всем ее атомам в последнем столбце.

Стоит отметитьчто файлы pdb имеют свой собственный формат, который необходимо поддерживать.Действительно, dest.pdb был создан строкой кода:

var=0.00
awk '{printf "%4s%7.0f%3s%6s%2s%4.0f%12.3f%8.3f%8.3f%6.2f%7.2f\n", 
$1, $2, $3, $4, $5, $6, $7, $8, $9, $10, $var}' < a.pdb >> dest.pdb

В моих неудачных попытках я пытался сопоставить целые числа в column1 source.csv с 6-м полем dest.pdb.Если это так, замените 11-е поле файла pdb значениями во втором столбце и в той же строке source.csv.Поскольку это блок кода гораздо большего скрипта bash, я попытался сделать это, используя только инструменты bash.В частности, я много пробовал с awk, что-то вроде:

while read -r n score; do
awk -v x=$n -v y=$score '{if ($6 == $x) printf "%4s%7.0f%3s%6s%2s%4.0f%12.3f%8.3f%8.3f%6.2f%7.2f\n", $1, $2, $3, $4, $5, $6, $7, $8, $9, $10, $y}' dest.pdb 
done < source.csv >> results.pdb

Я застрял на этом этапе и был бы признателен за любую помощь, чтобы двигаться дальше.Спасибо.

ОБНОВЛЕНИЕ. Я решил вышеуказанные проблемы с помощью следующей стратегии: (1) Удалите ненужные атомы водорода, которые запутались, потому что они имеют трехзначные имена атомов.Я сделал с помощью специальной небольшой программы, которая называется Reduce [http://kinemage.biochem.duke.edu/software/reduce.php][1] (2), удалите строки заголовков, вставленные в Reduce:

awk '$1 ~ /^ATOM/' A_temp.pdb >> A.pdb

(3) используйте read как парсер в bash:

while  read -r atom anum aname resname chain resnum x y z q b ; do
       while  read -r n evol ; do
       if [[ ${resnum} == ${n} ]]
       then
       echo "$atom $anum $aname $resname $chain $resnum $x $y $z $q $evol" >> out.pdb
      fi
      done < source.tsv
     done < A.pdb

(4) Отформатируйте файл A.pdf с помощью awk:

awk '{printf "%4s%7.0f%5s%4s%2s%4.0f%12.3f%8.3f%8.3f%6.2f%6.2f\n", $1, $2, $3, $4, $5, $6, $7, $8, $9, $10, $11}' < out.pdb >> Aform.pdb

echo 'TER' >> Aform.pdb

Answer 1

Bash не является полностью выращенным идеальным инструментом для разбора.Bash хорошо работает с разделенными файлами и не очень хорош в этом.Я бы предложил написать простую программу на C ++ или python для вашей цели, она бы работала проще, быстрее и безопаснее.

Для bash нам нужно сначала преобразовать файл во что-нибудь пригодное для синтаксического анализа.Я решил использовать вкладку для разделения полей.Поля в формате , который вы разместили, имеют постоянную длину - поэтому я вставлю вкладки в указанные позиции в файле.«Указанные позиции» - это окончания / начала полей в строках файла.Я вставляю вкладки с конца, чтобы нумерация символов не путалась.

sed '
    s/./\t&/79;
    s/./\t&/77;
    s/./\t&/61;
    s/./\t&/55;
    s/./\t&/47;
    s/./\t&/39;
    s/./\t&/31;
    s/./\t&/27;
    s/./\t&/23;
    s/./\t&/22;
    s/./\t&/18;
    s/./\t&/17;
    s/./\t&/13;
    s/./\t&/7;
' 

После этого мы можем использовать все стандартные утилиты Unix для разбора файлов.Я использовал просто объединить со строкой пользовательского формата, чтобы объединить файлы и распечатать то, что я хочу.Затем я могу просто удалить разделитель, чтобы восстановить формат.

# join on field 7 - resSeq
join -17 -21 -t$'\t' \
    -o 1.1,1.2,1.3,1.4,1.5,1.6,1.7,1.8,1.9,1.10,1.11,1.12,2.2 \
    -- <(
        # convert input file to tab separated file
        <dest.pdb \
        sed '
            s/./\t&/79;
            s/./\t&/77;
            s/./\t&/61;
            s/./\t&/55;
            s/./\t&/47;
            s/./\t&/39;
            s/./\t&/31;
            s/./\t&/27;
            s/./\t&/23;
            s/./\t&/22;
            s/./\t&/18;
            s/./\t&/17;
            s/./\t&/13;
            s/./\t&/7;
        ' |
        # sort it for join
        sort -t$'\t' -k7
    ) <(
        # convert input file for tab separated file
        <source.csv \
        tr ',' '\t' |
        # well squeese spaces, we format below
        tr -d ' ' |
        # the resSeq is 7 field length
        # you need to explicitly add spaces in front of the number
        # decide on which format do you use
        xargs -n2 printf "% 4d\t%7.2f\n" |
        # sort for join
        sort -t$'\t' -k1
) |
# I noticed the file was sorted using the `serial` field
sort -t$'\t' -n -k2 |
# and finally - remove bash separator
tr -d '\t'

Я предположил, что файл отсортирован по второму столбцу, поэтому я отсортировал его на выходе.Я также не знаю, как интерпретировать последнюю длину поля Real(6.2) - я использовал модификатор printf "%7.2", чтобы напечатать его.

Live версия доступна на tutorialspoint

Answer 2

Если я правильно понимаю данные, которые вы хотите, MET = 11.25 и GLU = 16.49

Я бы просто выполнил две команды sed из командной строки.Мой седь ржавый, но примерно (т.е. я его не проверял), если вы идете от 0,00 -> 11,25 и 0,00 -> 16,49,

sed 's/^(.+MET.+)0\.00/\1 11\.25/' myfile.pdb > outfile
sed 's/^(.+GLU.+)0\.00/\1 16\.49/' myfile.pdb >> outfile

Просто переверните числа между первымнабор скобок '/ first /' для второго набора скобок '/ first / second /', если вы хотите изменить операцию, т. е. с 11.25 -> 0.00 и т. д.

sed пользователи будутусилить это (и, вероятно, отладить это!), создав его как один скрипт.Я буду делать это только под принуждением.

В предлагаемых здесь решениях используется регулярное выражение, и, на мой взгляд, оно лучше реализовано на других языках, но sed - это простой способ выучить регулярное выражение, потому что оно чисторегулярное выражение и не требует никаких других знаний в области программирования.

Answer 3

Вы можете делать то, что вам нужно, используя ассоциативный массив для хранения значений csv и затем циклически перебирая строки в dest.pdb, вызывая sed, чтобы выполнить замену на основе 26-го символа в строкев качестве индекса для массива, используемого в подстановке sed.

( note: ) это нарушает правило эффективности вызова sed (или любой другой утилиты) в теле вашегоЦикл, но это необходимо в некоторых случаях. Если у вас есть миллион строк для преобразования - иди за чашкой кофе)

Причина в том, что это сложная задача, которую непросто сделать за один раз awkcall ваше второе поле содержит пробел , например

...    1  N   ...

(хотя вы можете включать условные выражения на основе NF, и это то, что вы в конечном итоге захотите сделать, следующее былочто первым пришло в голову)

Реализация в коротком сценарии будет такой:

#!/bin/bash

declare -A arr  ## declare & fill associative array from source.csv
while IFS=$' ,\t\n' read n v; do 
    arr[$n]=$v
done <source.csv

## read/covert last value in line using 26th char as index for array
#  and simple sed substitution to replace. 
#  (use braced group to redirect output)
{
while IFS= read -r line; do 
    sed "s/0[.]00$/${arr[${line:25:1}]}/" <<< "$line"
done <dest.pdb
} > results.pdb

Пример использования / Вывод

С вашими файлами данныхкак вы представили выше,простое выполнение скрипта в каталоге, содержащем source.csv и dest.pdb, приведет к results.pdb следующим образом:

$ cat results.pdb
ATOM      1  N   MET A   1     105.382 119.360 102.631  1.00   11.25
ATOM      2 CA   MET A   1     105.155 118.751 103.942  1.00   11.25
ATOM      3 HA   MET A   1     104.645 119.496 104.551  1.00   11.25
ATOM      4 CB   MET A   1     104.212 117.542 103.804  1.00   11.25
ATOM      5HB1   MET A   1     104.120 117.057 104.775  1.00   11.25
ATOM      6HB2   MET A   1     104.631 116.826 103.095  1.00   11.25
ATOM      7 CG   MET A   1     102.801 117.937 103.353  1.00   11.25
ATOM      8HG1   MET A   1     102.862 118.327 102.336  1.00   11.25
ATOM      9HG2   MET A   1     102.436 118.736 103.999  1.00   11.25
ATOM     10 SD   MET A   1     101.579 116.590 103.371  1.00   11.25
ATOM     11 CE   MET A   1     101.404 116.275 105.156  1.00   11.25
ATOM     12HE1   MET A   1     100.603 115.555 105.325  1.00   11.25
ATOM     13HE2   MET A   1     102.327 115.865 105.565  1.00   11.25
ATOM     14HE3   MET A   1     101.158 117.201 105.676  1.00   11.25
ATOM     15  C   MET A   1     106.423 118.387 104.697  1.00   11.25
ATOM     16  O   MET A   1     107.511 118.334 104.134  1.00   11.25
ATOM     17  N   GLU A   2     106.296 118.095 105.999  1.00   16.49
ATOM     18  H   GLU A   2     105.398 118.148 106.454  1.00   16.49
ATOM     19 CA   GLU A   2     107.495 117.802 106.786  1.00   16.49
ATOM     20 HA   GLU A   2     108.068 118.718 106.664  1.00   16.49
ATOM     21 CB   GLU A   2     107.242 117.714 108.295  1.00   16.49
ATOM     22HB1   GLU A   2     106.839 116.732 108.520  1.00   16.49
ATOM     23HB2   GLU A   2     106.494 118.455 108.581  1.00   16.49
ATOM     24 CG   GLU A   2     108.519 117.970 109.128  1.00   16.49
ATOM     25HG1   GLU A   2     108.323 117.660 110.155  1.00   16.49
ATOM     26HG2   GLU A   2     109.328 117.336 108.762  1.00   16.49
ATOM     27 CD   GLU A   2     109.002 119.432 109.126  1.00   16.49
ATOM     28OE1   GLU A   2     109.449 119.916 108.058  1.00   16.49
ATOM     29OE2   GLU A   2     109.026 120.057 110.206  1.00   16.49
ATOM     30  C   GLU A   2     108.446 116.757 106.163  1.00   16.49
ATOM     31  O   GLU A   2     109.650 117.015 106.154  1.00   16.49

Просмотрите все и дайте мне знать, если у вас есть вопросы.Я пойду посмотрю на дополнительную логику с awk и дам вам знать, если я что-нибудь там придумаю.