根据ID从FASTA文件中批量提取序列【Python】

根据ID从FASTA文件中批量提取序列【Python】

生信问题记录

我的需求

input：

1. FASTA文件，含六千余个蛋白序列。命名为FA.fasta
2. txt文件，经过interpro注释后，筛选出五千余个蛋白，将五千余个蛋白ID导出到txt文件中，每行一个。命名为ID.txt
   
   

output：

• FASTA文件，根据ID.txt里的ID从FA.fasta里提取出对应的五千余个序列。提取结果也以文件的格式保存为out_file.fasta

环境

Windows 10 64x
命令行运行py

解决流程

如果编程娴熟，可以直接用perl或Python写脚本；碍于个人技艺生疏，直接面向百度编程。搜索发现，解决方案非常多，最终某位网友的Python脚本成功运行

提示 如果你是用ID行，形如 >RLS74643.1 dihydrofolate reductase [Planctomycetes bacterium] 来提取序列，那么更为简单，可以使用联川生物的免费在线小工具

参考

1c43f522e1c3的帖子

https://www.jianshu.com/p/b7032dfae44e

步骤

直接复制网友的代码会出很多错，要先安装click模块和biopython的包

1. 安装包

pip install click
pip install biopython

2. 复制粘贴，另存为get_seqs_by_id.py

#coding:utf-8

import click

from Bio import SeqIO

@click.command()

@click.option('-f', '--fastafile', help='Input a fasta file', required=True)

@click.option('-i', '--idfile', help='Input an idlist', required=True)

@click.option('-o', '--outfile', help='Input the name of result file', default='result_out.fa')

def main(fastafile="FA.fa",idfile="ID.txt",outfile= "result_out.fa"):

    with open(idfile) as id_handle:
        wanted = set(line.rstrip("\n").split(None,1)[0] for line in id_handle)

    print("Found %i unique identifiers in %s" % (len(wanted), idfile))

    records = (r for r in SeqIO.parse(fastafile, "fasta") if r.id in wanted)

    count = SeqIO.write(records, outfile, "fasta")

    print("Saved %i records from %s to %s" % (count, fastafile, outfile))

    if count < len(wanted):

        print("Warning %i IDs not found in %s" % (len(wanted) - count, fastafile))
 
if __name__ == '__main__':
    main()

3. 将输入文件和py文件放在同一目录下，命令行输入

python get_seqs_by_id.py -f **.fasta -i ID.txt -o out_file.fasta

如果顺利出结果再好不过了，然而我报了错

报错不要慌，只有最后一行是关键信息

debug

错误的意思是：Unicode的解码（Decode）出现错误（Error）了，以gbk编码的方式去解码（该字符串变成Unicode），但是此处通过gbk的方式，却无法解码（can’t decode ）。“illegal multibyte sequence”意思是非法的多字节序列，即没法（解码）了。

于是我在open函数这一行加入encoding=‘UTF-8’

with open(idfile,encoding='UTF-8') as id_handle:

结果又双叒报错了！

耐心搜索，发现可能是用Python读文件（txt或者csv），出现编码错误
打开我的ID.txt一看，从Excel复制的信息居然是UTF-16！怪不得一直报错。把txt文件用UTF-8另存就好了……

运行成功

运行成功是这样的，输出一个fasta文件，正是我想要的。