ゲノムファイル全体から1つの周波数行列を計算するには？

Question

私は単純に、このようなパターンを含むHUGEファイル内の単一ヌクレオチド頻度（A、T、C、G）を計算しようとしています：TTTGTATAAGAAAAAATAGG。私のようなファイル全体の出力の1行与えるだろう

：ここ

The single nucleotide frequency matrix of T.volcanium Genome is: {'A': [234235], 'C': [234290], 'G': [32456], 'T': [346875]} 

は私の（ファイルパスのない、オープン、クローズおよびメイン）コードである

---

def freq_dict_of_lists_v1(dna_list): 
    n = max([len(dna) for dna in dna_list]) 
    frequency_matrix = { 
     'A': [0] * n, 
     'C': [0] * n, 
     'G': [0] * n, 
     'T': [0] * n, 
    } 
    for dna in dna_list: 
     for index, base in enumerate(dna): 
      frequency_matrix[base][index] += 1 

    return frequency_matrix 

for line in file: 
    dna_list = file.readline().rstrip("\n") 
    frequency_matrix = freq_dict_of_lists_v1(dna_list) 
    print("The single nucleotide frequency matrix of T.volcanium Genome is: ") 
    pprint.pprint(frequency_matrix) 

---

これは私の出力です。

The single nucleotide frequency matrix of T.volcanium Genome is: 
{'A': [21], 'C': [10], 'G': [11], 'T': [18]} 
The single nucleotide frequency matrix of T.volcanium Genome is: 
{'A': [31], 'C': [6], 'G': [4], 'T': [19]} 
The single nucleotide frequency matrix of T.volcanium Genome is: 
{'A': [23], 'C': [9], 'G': [10], 'T': [18]} 
The single nucleotide frequency matrix of T.volcanium Genome is: 
{'A': [17], 'C': [8], 'G': [9], 'T': [26]} 
The single nucleotide frequency matrix of T.volcanium Genome is: 
{'A': [15], 'C': [13], 'G': [9], 'T': [23]} 
The single nucleotide frequency matrix of T.volcanium Genome is: 
{'A': [21], 'C': [12], 'G': [10], 'T': [17]} 
The single nucleotide frequency matrix of T.volcanium Genome is: 
{'A': [20], 'C': [9], 'G': [12], 'T': [19]} 
The single nucleotide frequency matrix of T.volcanium Genome is: 
{'A': [15], 'C': [15], 'G': [10], 'T': [20]} 
The single nucleotide frequency matrix of T.volcanium Genome is: 
{'A': [20], 'C': [11], 'G': [10], 'T': [19]} 
The single nucleotide frequency matrix of T.volcanium Genome is: 
{'A': [26], 'C': [13], 'G': [7], 'T': [14]} 
The single nucleotide frequency matrix of T.volcanium Genome is: 
{'A': [12], 'C': [13], 'G': [13], 'T': [22]} 
The single nucleotide frequency matrix of T.volcanium Genome is: 
{'A': [20], 'C': [16], 'G': [9], 'T': [15]} 
The single nucleotide frequency matrix of T.volcanium Genome is: 
{'A': [22], 'C': [12], 'G': [6], 'T': [20]} 

したがって、1行に1行ずつ計算しています。 私はforループを取り出すか、readlinesを離そうとしましたが、ファイル内の1行だけの出力が1行しか得られません。ファイル全体ではありません。

私はこれを考え直しているように感じます。私は、ファイル全体を読んで、1つの出力を合計周波数で出力する簡単な方法があると確信しています...どんな洞察も感謝しています。

Answer 1

解決策に2つの問題があります。

1. あなたの質問に、それはあなたがあなたが一度行あたりの関数を呼び出しているすべての行
2. 渡ってカウントを追跡したいと言うときは、位置ごとに拠点を追跡しています。

私の編集内容は次のとおりです。説明のためのコメントを参照してください

def freq_dict_of_lists_v1(dna_list): 
    frequency_matrix = { # We are only keeping one variable per base 
     'A': [0],   # so that we calculate counts across all lines 
     'C': [0], 
     'G': [0], 
     'T': [0], 
    } 
    for dna in dna_list: 
     for base in dna: # No longer need index, so I removed enumerate 
      frequency_matrix[base] += 1 # Change here since dict structure changed 

    return frequency_matrix 

# Unlike before, we are now appending all the lines into dna_list 
for line in file: 
    dna_list.append(file.readline().rstrip("\n")) 

# Calling freq_dict_of_lists_v1 on ALL the lines at once (it is now out of loop) 
frequency_matrix = freq_dict_of_lists_v1(dna_list) 
print("The single nucleotide frequency matrix of T.volcanium Genome is: ") 
pprint.pprint(frequency_matrix) 

この解決策の注意点は、ファイル内のすべての塩基が大文字であることを確認することです。また、ACGT以外の文字がないことを確認してください（シーケンスには特殊な空白文字などがあります）。他の文字がある場合はthis threadを参照してください。デフォルトのエントリはGapのようになります。

Answer 2

MBが何を意味するのかわかりませんか？ GB？が、これは最も簡単な解決策です。ただし、ファイル全体がメモリにロードされることに注意してください。

# open file with sequence 
with open(path_to_file) as f: 
    seq = f.read() 

# count element A in sequence 
seq.count('A')