0

私は、Python を少し学ぼうとする最初の一歩を踏み出したところです。現在、バイオインフォマティクスの Python スキルを教えることを目的とした Rosalind オンライン コースに取り組んでいます。(ちなみに非常に良いです。参照: rosalind.info)

私は1つの特定の問題に苦しんでいます。次の形式の FASTA 形式のファイルがあります。

>Sequence_Header_1
ACGTACGTACGTACGTACGT
ACGTACGTACGTACGTACGT
>Sequence_Header_2
ACGTACGTACGTACGTACGT
ACGTACGTACGTACGTACGT

ファイルの各エントリ (ヘッダーを除く) で G と C の割合を計算し、この数値を返す必要があります。例:

>Sequence_Header_1
48.75%
>Sequence_header_2
52.43%

これまでの私のコードは次のとおりです。

file = open("input.txt" , "r")
for line in file:
    if line.startswith(">"):
        print(line.rstrip())        
    else:
        print ('%3.2f' % (line.count('G')+line.count('C')/len(line)*100))
file.close()

これは、私が必要とするほとんどのことを行っています。シーケンスデータが複数行にまたがる場所で問題が発生しています。現時点では、各エントリに対して単一の数値を返すのではなく、ファイル内のすべての行の % GC コンテンツを取得しています。例:

>Sequence_Header_1
48.75%
52.65%
>Sequence_header_2
52.43%
50.25%

複数の行にまたがるデータに数式を適用するにはどうすればよいですか?

前もって感謝します、

4

5 に答える 5

1

あなたの質問に対する直接的な回答ではありませんが、より良い方法だと思います。Python でさらにバイオインフォマティクスを行う予定がある場合は、biopython をご覧ください。fasta ファイルやその他の一般的なシーケンス操作 (およびそれ以上の操作) を処理します。

たとえば、次のようになります。

from Bio import SeqIO
from Bio.SeqUtils import GC

for rec in SeqIO.parse("input.txt", "fasta"):
    print rec.id,GC(rec.seq)
于 2013-07-11T15:56:03.760 に答える
0

次のように、fasta 形式を解析して、>ID をキー、シーケンスを値として辞書を作成できます。

    from collections import defaultdict

    def parse_fasta(dataset):
        "Parses data in FASTA format, returning a dictionary of ID's and values"
        records = defaultdict(str)
        record_id = None
        for line in [l.strip() for l in dataset.splitlines()]:
            if line.startswith('>'):
                record_id = line[1:]
            else:
                records[record_id] += line
        return records

または、このコードを少し書き直してタプル/リストを作成することもできます。辞書は既に索引付けされているため、私は辞書を好みます。それでも助けが必要な場合は、ロザリンドのサイトで私を見つけることができます.

于 2014-02-22T01:40:57.823 に答える
0

これがあなたが探しているものだと思います：

# Read the input file
with open("input.txt", "r") as f:
    s = f.read()

# Loop through the Sequences
for i, b in enumerate(s.split("Sequence_Header_")):
    if not b: continue # At least the first one is empty 
                       # because there is no data before the first sequence
    # Join the lines with data
    data = "".join(b.split("\n")[1:])

    # Print the result
    print("Sequence_Header_{i}\n{p:.2f}%".format(
        i=i, p=(data.count('G')+data.count('C'))/len(data)*100))

注：あなたの例では「>」記号が見つかりません。ヘッダーが > で始まる場合は、コードを s.split(">") に書き直すことができ、コードは問題ないはずです。

于 2013-07-10T14:43:10.950 に答える
0

ファイル全体をメモリに保持できない場合があります。すべてを一度に行うことはできないと仮定するとs = f.read()、新しいシーケンスが開始されるまで、文字数と総文字数の実行中のカウントを保持する必要があります。このようなもの：

file = open("input.txt" , "r")
# for keeping count:
char_count = 0
char_total = 0
for line in file:
    if line.startswith(">"):
        if char_total != 0:
            # starting a new sequence, so calculate pct for last sequence
            char_pct = (char_count / char_total) * 100
            print ('%3.2f' % char_pct)
            # reset the count for the next sequence
            char_total = 0
            char_count = 0
        print(line.rstrip())        
    else:
        # continuing a sequence, so add to running counts
        char_count += line.count('G') + line.count('C')
        char_total += len(line)
file.close()
于 2013-07-10T14:48:49.577 に答える