将文本添加到 fasta 文件的标头

abr*_*ham 2 bash awk r fasta str-replace

我有一个基因组的 fasta 文件(txt),类似于:

$ cat Strain-01.faa
>IMEHDJCA_03186 Serine/threonine-protein phosphatase 2
MEFKHRFIDGSRYQRIFVIGDIHGKLALLQDTLKRVDFHGERDLLISVGDLIDRGPDSVG
VLDYYQTHDWFEAVMGNHEWMMVNALDAQNKLERSEKEAYFIKIWHRNGCEWSQNL
>IMEHDJCA_03187 Serine transporter
MKESRETLNFSDTLPTETWTKHDTHWVLSLFGTAVGAGILFLPINLGIGGFWPLVLLALL
AFPMTFWGHRALARFVLSSKQADADFTDVVEEHFGAKAGRLISLLYFLSIFPILLIYGVG
>IMEHDJCA_03189 hypothetical protein
MNNQRHGITFGIERIGSQTILVFKATGTLTHQDYQAIAPVLEAALAGINRQQMNMLADIS
EFSGWEPRAAWDDFQLGLKIGFSVNKVAVYGDKNWQELAAKVGSWFISGEMKSFGD
Run Code Online (Sandbox Code Playgroud)

我想添加一个基于 file.txt 中的列表的额外 ID。

$ cat file.txt
ID      Gene        Strain-01       Strain-02       Strain-03
ID_01   pphB        IMEHDJCA_03186  DIBHEKPI_01648  LLMDBGDK_00598
ID_02   group_1001  IMEHDJCA_03187  DIBHEKPI_01635  LLMDBGDK_00611
ID_03   group_1002  IMEHDJCA_03189  DIBHEKPI_01628  LLMDBGDK_00616
Run Code Online (Sandbox Code Playgroud)

例如,对于 fastaStrain-01.faa文件,其IMEHDJCA_03186id 对应于Strain-01,因此我想将ID_01列 ID 的编号(来自file.txt)添加到序列的标题中,如下所示:

  • ID_01对应于IMEHDJCA_03186
  • ID_02对应于IMEHDJCA_03187
  • ID_03对应于IMEHDJCA_03189

结果会是这样的:

$cat Strain-01_edited.faa
>ID_01 IMEHDJCA_03186 Serine/threonine-protein phosphatase 2
MEFKHRFIDGSRYQRIFVIGDIHGKLALLQDTLKRVDFHGERDLLISVGDLIDRGPDSVG
VLDYYQTHDWFEAVMGNHEWMMVNALDAQNKLERSEKEAYFIKIWHRNGCEWSQNL
>ID_02 IMEHDJCA_03187 Serine transporter
MKESRETLNFSDTLPTETWTKHDTHWVLSLFGTAVGAGILFLPINLGIGGFWPLVLLALL
AFPMTFWGHRALARFVLSSKQADADFTDVVEEHFGAKAGRLISLLYFLSIFPILLIYGVG
>ID_03 IMEHDJCA_03189 hypothetical protein
MNNQRHGITFGIERIGSQTILVFKATGTLTHQDYQAIAPVLEAALAGINRQQMNMLADIS
EFSGWEPRAAWDDFQLGLKIGFSVNKVAVYGDKNWQELAAKVGSWFISGEMKSFGD
Run Code Online (Sandbox Code Playgroud)

我只想将 ID 代码添加file.txt到 fasta 文件的标头中。

有想法这样做吗?或bashR或任何其他方式?

非常感谢

mar*_*rkp 6

更新#2:通过以下方式消除特定菌株(名称)处理awk

  • 我们将把所有可能的应变/ID 映射加载到awk
  • 这将允许处理任何*.faa文件,而无需知道菌株名称
  • 这将允许处理*.faa混合菌株的文件(不知道这是否是OP必须解决的问题)
  • 降低awk代码的复杂性(与UPDATE #1相比),但需要额外的内存来存储更多id[]数组条目

样本数据(第一个字段中的菌株混合):

# for this (nonsensical?) file the first 3 blocks include a strain
# from each of the 3 columns (of strain names) from file.txt; the
# 4th block contains a nonsensical strain that doesn't exist in
# file.txt (ie, 4th line should not see an insertion of a ID value)

$ cat Strain-mixed.faa
>IMEHDJCA_03186 Serine transporter
MKESRETLNFSDTLPTETWTKHDTHWVLSLFGTAVGAGILFLPINLGIGGFWPLVLLALL
AFPMTFWGHRALARFVLSSKQADADFTDVVEEHFGAKAGRLISLLYFLSIFPILLIYGVG
>DIBHEKPI_01635 Serine/threonine-protein phosphatase 2
MEFKHRFIDGSRYQRIFVIGDIHGKLALLQDTLKRVDFHGERDLLISVGDLIDRGPDSVG
VLDYYQTHDWFEAVMGNHEWMMVNALDAQNKLERSEKEAYFIKIWHRNGCEWSQNL
>LLMDBGDK_00616 hypothetical protein
MNNQRHGITFGIERIGSQTILVFKATGTLTHQDYQAIAPVLEAALAGINRQQMNMLADIS
EFSGWEPRAAWDDFQLGLKIGFSVNKVAVYGDKNWQELAAKVGSWFISGEMKSFGD
>NO_MATCH hypothetical protein
MNNQRHGITFGIERIGSQTILVFKATGTLTHQDYQAIAPVLEAALAGINRQQMNMLADIS
EFSGWEPRAAWDDFQLGLKIGFSVNKVAVYGDKNWQELAAKVGSWFISGEMKSFGD
Run Code Online (Sandbox Code Playgroud)

awk将所有菌株加载到数组中的新代码id[]

awk '
NR==1   { next }                           # skip 1st line of 1st file
FNR==NR { for (i=3; i<=NF; i++)            # for rest of 1st file load id[] with ...
              id[$i]=$1                    # all strain/ID combos
          next
        }
/^>/    {                                   # for 2nd file, if 1st column is ">"
          ndx=substr($1,2)                  # strip off ">"
          if ( ndx in id )                  # if 1st field (sans ">") is an index in id[] then ...
             ( $1=">" id[ndx] " " ndx )     # rewrite 1st field to include our id[] value
        }
1                                           # print current line (of 2nd file)
' file.txt Strain-mixed.faa
Run Code Online (Sandbox Code Playgroud)

这会生成:

>ID_01 IMEHDJCA_03186 Serine transporter
MKESRETLNFSDTLPTETWTKHDTHWVLSLFGTAVGAGILFLPINLGIGGFWPLVLLALL
AFPMTFWGHRALARFVLSSKQADADFTDVVEEHFGAKAGRLISLLYFLSIFPILLIYGVG
>ID_02 DIBHEKPI_01635 Serine/threonine-protein phosphatase 2
MEFKHRFIDGSRYQRIFVIGDIHGKLALLQDTLKRVDFHGERDLLISVGDLIDRGPDSVG
VLDYYQTHDWFEAVMGNHEWMMVNALDAQNKLERSEKEAYFIKIWHRNGCEWSQNL
>ID_03 LLMDBGDK_00616 hypothetical protein
MNNQRHGITFGIERIGSQTILVFKATGTLTHQDYQAIAPVLEAALAGINRQQMNMLADIS
EFSGWEPRAAWDDFQLGLKIGFSVNKVAVYGDKNWQELAAKVGSWFISGEMKSFGD
>NO_MATCH hypothetical protein
MNNQRHGITFGIERIGSQTILVFKATGTLTHQDYQAIAPVLEAALAGINRQQMNMLADIS
EFSGWEPRAAWDDFQLGLKIGFSVNKVAVYGDKNWQELAAKVGSWFISGEMKSFGD
Run Code Online (Sandbox Code Playgroud)

注意:此最新更新将为文件中的所有行执行 ID 插入Strain-{01,02}.faa(请参阅下面的更新 #1)。



更新#1:扩展原始答案以解决(我认为)Paul Hodges 关于概括答案以支持可变菌株名称的问题:

  • 动态确定使用哪一列菌株file.txt
  • 动态处理匹配的<strain>.faa文件

样本数据:

$ cat Strain-01.faa
>IMEHDJCA_03186 Serine/threonine-protein phosphatase 2
MEFKHRFIDGSRYQRIFVIGDIHGKLALLQDTLKRVDFHGERDLLISVGDLIDRGPDSVG
VLDYYQTHDWFEAVMGNHEWMMVNALDAQNKLERSEKEAYFIKIWHRNGCEWSQNL
>IMEHDJCA_03187 Serine transporter
MKESRETLNFSDTLPTETWTKHDTHWVLSLFGTAVGAGILFLPINLGIGGFWPLVLLALL
AFPMTFWGHRALARFVLSSKQADADFTDVVEEHFGAKAGRLISLLYFLSIFPILLIYGVG
>IMEHDJCA_03189 hypothetical protein
MNNQRHGITFGIERIGSQTILVFKATGTLTHQDYQAIAPVLEAALAGINRQQMNMLADIS
EFSGWEPRAAWDDFQLGLKIGFSVNKVAVYGDKNWQELAAKVGSWFISGEMKSFGD

# for this next file I simply copied data from OP's Strain-01.faa and
# modified the initial field for blocks 1 & 3; net result is we should
# see 2 of the blocks receive insertions of ID values

$ cat Strain-02.faa
>DIBHEKPI_01635 Serine/threonine-protein phosphatase 2
MEFKHRFIDGSRYQRIFVIGDIHGKLALLQDTLKRVDFHGERDLLISVGDLIDRGPDSVG
VLDYYQTHDWFEAVMGNHEWMMVNALDAQNKLERSEKEAYFIKIWHRNGCEWSQNL
>IMEHDJCA_03187 Serine transporter
MKESRETLNFSDTLPTETWTKHDTHWVLSLFGTAVGAGILFLPINLGIGGFWPLVLLALL
AFPMTFWGHRALARFVLSSKQADADFTDVVEEHFGAKAGRLISLLYFLSIFPILLIYGVG
>DIBHEKPI_01648 hypothetical protein
MNNQRHGITFGIERIGSQTILVFKATGTLTHQDYQAIAPVLEAALAGINRQQMNMLADIS
EFSGWEPRAAWDDFQLGLKIGFSVNKVAVYGDKNWQELAAKVGSWFISGEMKSFGD

$ cat Strain-XX.faa
cat: Strain-XX.faa: No such file or directory
Run Code Online (Sandbox Code Playgroud)

对原始答案进行一些修改awk,并包装在 ( bash)for循环中以处理不同的菌株:

for strain in Strain-01 Strain-02 Strain-XX
do
    printf "\n############### ${strain} / ${strain}.faa\n\n"

    awk -v strain="${strain}" '                 # pass bash variable in as awk variable (same name)

    NR==1   { for (i=3; i<=NF; i++)             # 1st row of 1st file, look for matching strain name
                  { if ( $i == strain )         # if we find a match then ...
                       { strain_ndx=i           # make note of the column and ...
                         next                   # skip to next line from 1st file
                       }
                  }

              # if we got here we did not find a matching strain name so 
              # print a message and exit from our awk script 

              print "Unable to locate entry for "strain" in "FILENAME". Aborting."
              exit
            }

    FNR==NR { id[$(strain_ndx)]=$1              # for rest of 1st file build array of ids
              next
            }

    /^>/    {                                   # for 2nd file, if 1st column is ">"
              ndx=substr($1,2)                  # strip off ">"
              if ( ndx in id )                  # if 1st field (sans ">") is an index in id[] then ...
                 ( $1=">" id[ndx] " " ndx )     # rewrite 1st field to include our id[] value
            }
    1                                           # print current line (of 2nd file)
    ' file.txt "${strain}.faa"
done
Run Code Online (Sandbox Code Playgroud)

这会生成:

############### Strain-01 / Strain-01.faa

>ID_01 IMEHDJCA_03186 Serine/threonine-protein phosphatase 2
MEFKHRFIDGSRYQRIFVIGDIHGKLALLQDTLKRVDFHGERDLLISVGDLIDRGPDSVG
VLDYYQTHDWFEAVMGNHEWMMVNALDAQNKLERSEKEAYFIKIWHRNGCEWSQNL
>ID_02 IMEHDJCA_03187 Serine transporter
MKESRETLNFSDTLPTETWTKHDTHWVLSLFGTAVGAGILFLPINLGIGGFWPLVLLALL
AFPMTFWGHRALARFVLSSKQADADFTDVVEEHFGAKAGRLISLLYFLSIFPILLIYGVG
>ID_03 IMEHDJCA_03189 hypothetical protein
MNNQRHGITFGIERIGSQTILVFKATGTLTHQDYQAIAPVLEAALAGINRQQMNMLADIS
EFSGWEPRAAWDDFQLGLKIGFSVNKVAVYGDKNWQELAAKVGSWFISGEMKSFGD

############### Strain-02 / Strain-02.faa

>ID_02 DIBHEKPI_01635 Serine/threonine-protein phosphatase 2
MEFKHRFIDGSRYQRIFVIGDIHGKLALLQDTLKRVDFHGERDLLISVGDLIDRGPDSVG
VLDYYQTHDWFEAVMGNHEWMMVNALDAQNKLERSEKEAYFIKIWHRNGCEWSQNL
>IMEHDJCA_03187 Serine transporter
MKESRETLNFSDTLPTETWTKHDTHWVLSLFGTAVGAGILFLPINLGIGGFWPLVLLALL
AFPMTFWGHRALARFVLSSKQADADFTDVVEEHFGAKAGRLISLLYFLSIFPILLIYGVG
>ID_01 DIBHEKPI_01648 hypothetical protein
MNNQRHGITFGIERIGSQTILVFKATGTLTHQDYQAIAPVLEAALAGINRQQMNMLADIS
EFSGWEPRAAWDDFQLGLKIGFSVNKVAVYGDKNWQELAAKVGSWFISGEMKSFGD

############### Strain-XX / Strain-XX.faa

Unable to locate entry for Strain-XX in file.txt. Aborting.
Run Code Online (Sandbox Code Playgroud)

原答案

一个awk想法:

awk '
FNR==NR { id[$3]=$1 ; next }                # for 1st file build array of ids
/^>/    {                                   # for 2nd file, if 1st column is ">"
          ndx=substr($1,2)                  # strip off ">"
          if ( ndx in id )                  # if 1st field (sans ">") is an index in id[] then ...
             ( $1=">" id[ndx] " " ndx )     # rewrite 1st field to include our id[] value
        }
1                                           # print current line (of 2nd file)
' file.txt fasta.dat
Run Code Online (Sandbox Code Playgroud)

对于给定的样本数据,这会生成:

>ID_01 IMEHDJCA_03186 Serine/threonine-protein phosphatase 2
MEFKHRFIDGSRYQRIFVIGDIHGKLALLQDTLKRVDFHGERDLLISVGDLIDRGPDSVG
VLDYYQTHDWFEAVMGNHEWMMVNALDAQNKLERSEKEAYFIKIWHRNGCEWSQNL
>ID_02 IMEHDJCA_03187 Serine transporter
MKESRETLNFSDTLPTETWTKHDTHWVLSLFGTAVGAGILFLPINLGIGGFWPLVLLALL
AFPMTFWGHRALARFVLSSKQADADFTDVVEEHFGAKAGRLISLLYFLSIFPILLIYGVG
>ID_03 IMEHDJCA_03189 hypothetical protein
MNNQRHGITFGIERIGSQTILVFKATGTLTHQDYQAIAPVLEAALAGINRQQMNMLADIS
EFSGWEPRAAWDDFQLGLKIGFSVNKVAVYGDKNWQELAAKVGSWFISGEMKSFGD
Run Code Online (Sandbox Code Playgroud)