【问题标题】:How do I convert the three letter amino acid codes to one letter code with python or R?如何使用 python 或 R 将三个字母的氨基酸代码转换为一个字母的代码?
【发布时间】:2012-10-06 13:39:20
【问题描述】:

我有一个如下所示的 fasta 文件。我想将three letter codes 转换为一个字母代码。我怎样才能用 python 或 R 做到这一点?

>2ppo
ARGHISLEULEULYS
>3oot
METHISARGARGMET

想要的输出

>2ppo
RHLLK
>3oot
MHRRM

您的建议将不胜感激!

【问题讨论】:

  • ARGHISLEULEULYS如何转换为RHLLK?逻辑是什么?
  • @Tichodroma:ARG = R,HIS = H,LEU = L,等等
  • @Junuxx 等?将完整的翻译列表添加到问题或至少链接到它会很有用。我想帮助解决这个问题,但除非我获得所有必要的信息,否则我无法提供帮助。
  • 我很好奇你在哪里找到这样一个文件 - 我从来没有见过使用这样的三个字母氨基酸代码的 FASTA 文件。

标签: python r bioinformatics biopython


【解决方案1】:

BioPython 已经内置了字典来帮助进行此类翻译。以下命令将显示可用字典的完整列表:

import Bio
help(Bio.SeqUtils.IUPACData)

您要查找的预定义字典:

Bio.SeqUtils.IUPACData.protein_letters_3to1['Ala']

【讨论】:

  • 这应该是选择的答案。一个小提示:在 Python3 中,至少该方法实际上是在模块 Bio.Data 下,而 Bio.SeqUtilis 从那里导入它,因此如果只想要当前命名空间中的方法 protein_letters_3to1,可以这样做:from Bio.Data.IUPACData import protein_letters_3to1
【解决方案2】:

用字典查单字代码:

d = {'CYS': 'C', 'ASP': 'D', 'SER': 'S', 'GLN': 'Q', 'LYS': 'K',
     'ILE': 'I', 'PRO': 'P', 'THR': 'T', 'PHE': 'F', 'ASN': 'N', 
     'GLY': 'G', 'HIS': 'H', 'LEU': 'L', 'ARG': 'R', 'TRP': 'W', 
     'ALA': 'A', 'VAL':'V', 'GLU': 'E', 'TYR': 'Y', 'MET': 'M'}

还有一个简单的函数来匹配整个字符串的三个字母代码和一个字母代码:

def shorten(x):
    if len(x) % 3 != 0: 
        raise ValueError('Input length should be a multiple of three')

    y = ''
    for i in range(len(x) // 3):
            y += d[x[3 * i : 3 * i + 3]]
    return y

测试你的例子:

>>> shorten('ARGHISLEULEULYS')
'RHLLK'

【讨论】:

  • 非常感谢您的回答。我是 python 新手。如何将输入文件解析为您的代码?
  • @user1725152:这取决于输入文件的格式。但我想它可能类似于for line in inputfile: print(shorten(line))
  • len(x) / 3 返回一个浮点数,因此如果您收到错误 TypeError: 'float' object cannot be interpreted as an integer 只需将其更改为: ``` for i in range(int(len(x)/3)): ```
  • @universvm:感谢您的评论。这是从 2012 年开始的,所以它是用 Python 2 编写的,len(x) / 3 将返回一个 int。更新了使用整数除法的答案。
【解决方案3】:

这是在 R 中执行此操作的一种方法:

# Variables:
foo <- c("ARGHISLEULEULYS","METHISARGARGMET")

# Code maps:
code3 <- c("Ala", "Arg", "Asn", "Asp", "Cys", "Glu", "Gln", "Gly", "His", 
"Ile", "Leu", "Lys", "Met", "Phe", "Pro", "Ser", "Thr", "Trp", 
"Tyr", "Val")
code1 <- c("A", "R", "N", "D", "C", "E", "Q", "G", "H", "I", "L", "K", 
"M", "F", "P", "S", "T", "W", "Y", "V")

# For each code replace 3letter code by 1letter code:
for (i in 1:length(code3))
{
    foo <- gsub(code3[i],code1[i],foo,ignore.case=TRUE)
}

结果:

> foo
[1] "RHLLK" "MHRRM"

请注意,我更改了变量名,因为变量名在 R 中不允许以数字开头。

【讨论】:

  • 这不好。以 TRPHISGLU 为例,您希望算法翻译如下 {TRP}{HIS}{GLU} -> WHE 但您的算法真正发生的是 TRP{HIS}{GLU} -> TR{PHE} -> TRF .您确实需要将foo 拆分为三个字符的子字符串以避免这种可能的交互。
【解决方案4】:
>>> src = "ARGHISLEULEULYS"
>>> trans = {'ARG':'R', 'HIS':'H', 'LEU':'L', 'LYS':'K'}
>>> "".join(trans[src[x:x+3]] for x in range(0, len(src), 3))
'RHLLK'

您只需将其余条目添加到trans dict。

编辑:

要制作trans 的其余部分,您可以这样做。档案table

Ala A
Arg R
Asn N
Asp D
Cys C
Glu E
Gln Q
Gly G
His H
Ile I
Leu L
Lys K
Met M
Phe F
Pro P
Ser S
Thr T
Trp W
Tyr Y
Val V

阅读:

trans = dict((l.upper(), s) for l, s in
             [row.strip().split() for row in open("table").readlines()])

【讨论】:

    【解决方案5】:

    您可以尝试查看并安装 Biopython,因为您正在解析 .fasta 文件,然后转换为单字母代码。不幸的是,Biopython 只有函数 seq3(在包 Bio::SeqUtils 中),它与你想要的相反。 IDLE 中的示例输出:

    >>>seq3("MAIVMGRWKGAR*")
    >>>'MetAlaIleValMetGlyArgTrpLysGlyAlaArgTer'
    

    不幸的是,没有“seq1”功能(还...),但我认为这可能对您将来有所帮助。至于你的问题,Junuxx 是正确的。创建一个字典并使用 for 循环以三个为一组读取字符串并进行翻译。这是一个与他提供的功能类似的功能,它包罗万象,也可以处理小写。

    def AAcode_3_to_1(seq):
        '''Turn a three letter protein into a one letter protein.
    
        The 3 letter code can be upper, lower, or any mix of cases
        The seq input length should be a factor of 3 or else results
        in an error
    
        >>>AAcode_3_to_1('METHISARGARGMET')
        >>>'MHRRM'
    
        '''
        d = {'CYS': 'C', 'ASP': 'D', 'SER': 'S', 'GLN': 'Q', 'LYS': 'K',
         'ILE': 'I', 'PRO': 'P', 'THR': 'T', 'PHE': 'F', 'ASN': 'N', 
         'GLY': 'G', 'HIS': 'H', 'LEU': 'L', 'ARG': 'R', 'TRP': 'W', 'TER':'*',
         'ALA': 'A', 'VAL':'V', 'GLU': 'E', 'TYR': 'Y', 'MET': 'M','XAA':'X'}
    
        if len(seq) %3 == 0:
            upper_seq= seq.upper()
            single_seq=''
            for i in range(len(upper_seq)/3):
                single_seq += d[upper_seq[3*i:3*i+3]]
            return single_seq
        else:
            print("ERROR: Sequence was not a factor of 3 in length!")
    

    【讨论】:

    • 从下一个版本 Biopython 1.61 开始,您将能够使用 Bio.SeqUtils.seq1(如果您喜欢处于领先地位,也可以从 github 存储库运行)。
    【解决方案6】:

    Biopython 有一个很好的解决方案

    >>> from Bio.PDB.Polypeptide import *
    >>> three_to_one('ALA')
    'A'
    

    对于你的例子,我将通过这一行来解决它

    >>> from Bio.PDB.Polypeptide import *
    >>> str3aa = 'ARGHISLEULEULYS'
    >>> "".join([three_to_one(aa3) for aa3 in [ "".join(g) for g in zip(*(iter(str3aa),) * 3)]])
    >>> 'RHLLK'
    

    他们可能会批评我这种类型的一个班轮:),但在我内心深处,我仍然爱着 PERL。

    【讨论】:

      【解决方案7】:

      使用 R:

      convert <- function(l) {
      
        map <- c("A", "R", "N", "D", "C", "E", "Q", "G", "H", "I",
                 "L", "K", "M", "F", "P", "S", "T", "W", "Y", "V")
      
        names(map) <- c("ALA", "ARG", "ASN", "ASP", "CYS", "GLU", "GLN",
                        "GLY", "HIS", "ILE", "LEU", "LYS", "MET", "PHE",
                        "PRO", "SER", "THR", "TRP", "TYR", "VAL")
      
        sapply(strsplit(l, "(?<=[A-Z]{3})", perl = TRUE),
               function(x) paste(map[x], collapse = ""))
      }
      
      convert(c("ARGHISLEULEULYS", "METHISARGARGMET"))
      # [1] "RHLLK" "MHRRM"
      

      【讨论】:

      • +1 用于将字符串拆分为 3 个字符的子字符串的巧妙方法。它展示了一些关于正则表达式匹配如何工作的有趣内容。
      • @fodel 非常感谢您的回答。我有超过 1000 个序列。它在一个文本文件中。首先,我必须将此文件导入到 r 中,并且必须将三个字母代码更改为一个字母。我已经显示了所需的输出。如果可以,请帮助我。
      • 我向您展示的函数将序列向量作为输入。如何将 FASTA 文件读入 R 中的序列向量是另一个问题。快速谷歌搜索,我可以为你指出至少三个不同的包:Biostrings (readFASTA)seqinr (read.fasta)bio3d (read.fasta)
      【解决方案8】:

      另一种方法是使用 R 中的 seqinriPAC 包。

      # install.packages("seqinr")
      # source("https://bioconductor.org/biocLite.R")
      # biocLite("iPAC")
      
      library(seqinr)
      library(iPAC)
      
      #read in file
      fasta = read.fasta(file = "test_fasta.fasta", seqtype = "AA", as.string = T, set.attributes = F)
      #split string
      n = 3
      fasta1 = lapply(fasta,  substring(x,seq(1,nchar(x),n),seq(n,nchar(x),n)))
      #convert the three letter code for each element in the list 
      fasta2 = lapply(fasta1, function(x) paste(sapply(x, get.SingleLetterCode), collapse = ""))
      
      # > fasta2
      # $`2ppo`
      # [1] "RHLLK"
      #
      # $`3oot`
      # [1] "MHRRM"
      

      【讨论】:

        【解决方案9】:
        my %aa_hash=(
          Ala=>'A',
          Arg=>'R',
          Asn=>'N',
          Asp=>'D',
          Cys=>'C',
          Glu=>'E',
          Gln=>'Q',
          Gly=>'G',
          His=>'H',
          Ile=>'I',
          Leu=>'L',
          Lys=>'K',
          Met=>'M',
          Phe=>'F',
          Pro=>'P',
          Ser=>'S',
          Thr=>'T',
          Trp=>'W',
          Tyr=>'Y',
          Val=>'V',
          Sec=>'U',                       #http://www.uniprot.org/manual/non_std;Selenocysteine (Sec) and pyrrolysine (Pyl)
          Pyl=>'O',
        );
        
        
            while(<>){
                    chomp;
                    my $aa=$_;
                    warn "ERROR!! $aa invalid or not found in hash\n" if !$aa_hash{$aa};
                    print "$aa\t$aa_hash{$aa}\n";
            }
        

        使用此 perl 脚本将三元组 a.a 代码转换为单字母代码。

        【讨论】:

          【解决方案10】:

          对于那些在 2017 年及以后登陆这里的人:

          这是一个单行 Linux bash 命令,用于将蛋白质氨基酸三字母代码转换为文本文件中的单字母代码。我知道这不是很优雅,但我希望这对搜索相同内容并希望使用单行命令的人有所帮助。

          sed 's/ALA/A/g;s/CYS/C/g;s/ASP/D/g;s/GLU/E/g;s/PHE/F/g;s/GLY/G/g;s/HIS/H/g;s/HID/H/g;s/HIE/H/g;s/ILE/I/g;s/LYS/K/g;s/LEU/L/g;s/MET/M/g;s/ASN/N/g;s/PRO/P/g;s/GLN/Q/g;s/ARG/R/g;s/SER/S/g;s/THR/T/g;s/VAL/V/g;s/TRP/W/g;s/TYR/Y/g;s/MSE/X/g' < input_file_three_letter_code.txt > output_file_single_letter_code.txt
          

          上述原始问题的解决方案,作为单个命令行:

          sed 's/.\{3\}/& /g' | sed 's/ALA/A/g;s/CYS/C/g;s/ASP/D/g;s/GLU/E/g;s/PHE/F/g;s/GLY/G/g;s/HIS/H/g;s/HID/H/g;s/HIE/H/g;s/ILE/I/g;s/LYS/K/g;s/LEU/L/g;s/MET/M/g;s/ASN/N/g;s/PRO/P/g;s/GLN/Q/g;s/ARG/R/g;s/SER/S/g;s/THR/T/g;s/VAL/V/g;s/TRP/W/g;s/TYR/Y/g;s/MSE/X/g' | sed 's/ //g' < input_file_three_letter_code.txt > output_file_single_letter_code.txt
          

          解释:

          [1] sed 's/.\{3\}/&amp; /g' 将拆分序列。它会在每 3 个字母后添加一个空格。

          [2] 管道中的第二个 'sed' 命令将获取上面的输出并转换为单字母代码。添加任何非标准残基为s/XYZ/X/g; 到这个命令。

          [3] 第三个“sed”命令,sed 's/ //g' 将删除空格。

          【讨论】:

            【解决方案11】:

            Python 3 解决方案。

            在我的工作中,令人烦恼的部分是氨基酸代码可以参考PDB / mmCIF文件中经常出现的修改后的代码,例如

            'Tih'-->'A'.

            所以映射可以超过 22 对。 Python 中的第 3 方工具,如

            Bio.SeqUtils.IUPACData.protein_letters_3to1

            无法处理。我最简单的解决方案是使用http://www.ebi.ac.uk/pdbe-srv/pdbechem 查找映射并将异常映射添加到我自己的函数中的字典中,每当我遇到它们时。

            def three_to_one(three_letter_code):
                mapping = {'Aba':'A','Ace':'X','Acr':'X','Ala':'A','Aly':'K','Arg':'R','Asn':'N','Asp':'D','Cas':'C',
                       'Ccs':'C','Cme':'C','Csd':'C','Cso':'C','Csx':'C','Cys':'C','Dal':'A','Dbb':'T','Dbu':'T',
                       'Dha':'S','Gln':'Q','Glu':'E','Gly':'G','Glz':'G','His':'H','Hse':'S','Ile':'I','Leu':'L',
                       'Llp':'K','Lys':'K','Men':'N','Met':'M','Mly':'K','Mse':'M','Nh2':'X','Nle':'L','Ocs':'C',
                       'Pca':'E','Phe':'F','Pro':'P','Ptr':'Y','Sep':'S','Ser':'S','Thr':'T','Tih':'A','Tpo':'T',
                       'Trp':'W','Tyr':'Y','Unk':'X','Val':'V','Ycm':'C','Sec':'U','Pyl':'O'} # you can add more
                return mapping[three_letter_code[0].upper() + three_letter_code[1:].lower()]
            

            另一种解决方案是在线检索映射(但url和html模式可能会随着时间而改变):

            import re
            import urllib.request
            
            def three_to_one_online(three_letter_code):
                url = "http://www.ebi.ac.uk/pdbe-srv/pdbechem/chemicalCompound/show/" + three_letter_code
                with urllib.request.urlopen(url) as response:
                    single_letter_code = re.search('\s*<td\s*>\s*<h3>One-letter code.*</h3>\s*</td>\s*<td>\s*([A-Z])\s*</td>', response.read().decode('utf-8')).group(1)
                return single_letter_code
            

            为了简单起见,这里我直接使用 re 而不是 html 解析器。

            希望这些能有所帮助。

            【讨论】:

              猜你喜欢
              • 1970-01-01
              • 1970-01-01
              • 1970-01-01
              • 1970-01-01
              • 1970-01-01
              • 1970-01-01
              • 2022-01-23
              • 2017-09-29
              • 2021-04-23
              相关资源
              最近更新 更多