BiopythonでIDSに基づくFASTAファイルをフィルタリングする

Question

私はPythonでプログラミングするのがとても新しいです。私はいくつかの植物種のタンパク質配列を含むファーストファイルを持っています。

各配列に含まれるアミノ酸の数に基づいてそれらをフィルタリングしたいと思います。基準は20アミノ酸を超える配列である。

biopython cookbookのリソースで20を超えるアミノ酸配列を得ることができます。しかし、ファイルに書き込もうとすると、これは私にこのErrorを与えます。私はこのエラーを解決することができません。また、出力ファイルに各シーケンスのIDを持たせたいと思います。私を助けてください！

コード：

import Bio 
from Bio import SeqIO 
for s_record in SeqIO.parse('arabidopsis_thaliana_proteome.ath.tfa','fasta'): 
    name = s_record.id 
    seq = s_record.seq 
    seqLen = len(s_record) 
    if seqLen >20: 
     desired_proteins=seq 
     output_file=SeqIO.write(desired_proteins, "filtered.fasta","fasta") 
output_file 

入力ファイル：Arabidopsis Thaliana

>AT5G16970 
MTATNKQVILKDYVSGFPTESDFDFTTTTVELRVPEGTNSVLVKNLYLSCDPYMRIRMGKPDPSTAALAQAYTPGQPIQGYGVSRIIESGHPDYKKGDLLWGIVAWEEYSVITPMTHAHFKIQHTDVPLSYYTGLLGMPGMTAYAGFYEVCSPKEGETVYVSAASGAVGQLVGQLAKMMGCYVVGSAGSKEKVDLLKTKFGFDDAFNYKEESDLTAALKRCFPNGIDIYFENVGGKMLDAVLVNMNMHGRIAVCGMISQYNLENQEGVHNLSNIIYKRIRIQGFVVSDFYDKYSKFLEFVLPHIREGKITYVEDVADGLEKAPEALVGLFHGKNVGKQVVVVARE* 

>AT4G32100 
MATNACKFLCLVLLFAFVTQGYGDDSYSLESLSVIQSKTGNMVENKPEWEVKVLNSSPCYFTHTTLSCVRFKSVTPIDSKVLSKSGDTCLLGNGDSIHDISFKYVWDTSFDLKVVDGYIACS* 

は、事前にありがとう:)

Answer 1

をここBioPythonチュートリアルによると：

http://biopython.org/wiki/SeqIO

SeqIO.parse()への最初のパラメータはファイルハンドルではなくファイル名でなければなりません：

from Bio import SeqIO 
with open("example.fasta", "rU") as handle: 
    for record in SeqIO.parse(handle, "fasta"): 
     print(record.id) 

これは動作するはずです：

import Bio 
from Bio import SeqIO 
fh=open('arabidopsis_thaliana_proteome.ath.tfa') 
for s_record in SeqIO.parse(fh,'fasta'): 
    name = s_record.id 
    seq = s_record.seq 
    seqLen = len(s_record) 
    if seqLen >20: 
     desired_proteins=seq 
     output_file=SeqIO.write(desired_proteins, "filtered.fasta","fasta") 
output_file 
fh.close() 

Answer 2

bioawkを使用して溶液（バイオインフォマティクスのための便利なawkの修正版）代わりにbiopythonの：

bioawk -c fastx 'length($seq) > 20 {print ">"$name"\n"$seq}' arabidopsis_thaliana_proteome.ath.tfa > filtered.fasta 

-c fastx fastx/FASTQとしてファイルを解析するbioawkを伝えますフォーマット。これは、普通の'condition {action}' awk構文を使用して使用できるnameとseqの変数を定義します。

私はbioawkがほとんどのLinuxディストリビューション用のパッケージが、ソースからのインストールはあまりにも複雑ではありませんわからない： https://silico-sciences.com/2015/12/13/install-bioawk-on-ubuntu/

Answer 3

例（長さ300の下のシーケンスのフィルタリング）ほぼ正確にこれを行うことがあります、 "入力/出力の例 - シーケンス長によるフィルタリング" on http://biopython.org/wiki/SeqIO

これは、メモリー効率の高いジェネレータ表現に対して、メモリ内のレコードのリストを構築するオプションになります。

第三の選択肢は、出力ハンドルを使用してSeqIO.writeへの複数の呼び出しを行うことである。

from Bio import SeqIO 

with open("short_seqs.fasta", "w") as out_handle: 
    for record in SeqIO.parse("cor6_6.gb", "genbank"): 
     if len(record.seq) < 300: 
      SeqIO.write(record, out_handle, "fasta") 

しかし、このアプローチは、あなたが保つことができるFASTA、GenBankの、EMBL、のような単純なファイル形式にのみ適しています既存のファイルにレコードを追加します。これはXML出力フォーマットのようなものでは機能しません。リンクされたページに表示されるアプローチはより良いです（レコードのジェネレータ/イテレータを使用してSeqIO.writeを1回呼び出す）。