2014-12-02 94 views
0

我有一個程序,正在讀取文件「側翼seqs」包含文本列各不同的東西的意思:程序沒有給我輸出我想

1 1 44457990 TAA CTCTCCTAAAGGACC 
1 1 44461833 TGA CCAGCCTGAAGGGCT 
1 1 148594641 TAA CCACAATAAGCAGCT 
1 1 43241066 TGA ACTCACTGAGAGTGG 
1 1 43240880 TAG CTTCTCTAGGAATGG ... 

首先山坳:染色體數目,第二col:DNA鏈,第三列:DNA中終止密碼子的位置,第四列:終止密碼子,第五列:終止密碼子周圍的上游和下游6個鹼基,即每個終止密碼子的側翼序列。

現在,我的程序應該讀取此文件,並從側翼序列列中提取每個終止密碼子前後的3個鹼基,並寫入包含兩列的文件:終止密碼子,然後是側翼序列。該文件應該包含所有三個終止密碼子TAA,TAG和TGA的側翼序列,但是當我運行該程序時,它只給出了TGA終止密碼子的側翼序列,但是對於其餘的則沒有。

這裏是什麼樣的OUTFILE看起來像一個例子:

TGA GGGCTT 1 
TGA GAACGT 2 
TGA CTTCTT 17 
TGA CACCCT 15 
TGA GAACGG 1 
TGA GAACGC 3 

我看不到我要去的地方錯了,但我不是很有經驗,所以我敢肯定,我失去了一些東西簡單。我很感激任何幫助發現我的錯誤!下面是代碼:

bases = ['A','T','C','G'] 
sequenceCount = {} 
for x1 in bases: 
    for x2 in bases: 
     for x3 in bases: 
      for x4 in bases: 
       for x5 in bases: 
        for x6 in bases: 
         sequenceCount[x1+x2+x3+x4+x5+x6] = 0 
infile = open('flanking seqs.txt','rU') 
outfile = open('context resultsNEW.txt','w') 

for line in infile: 
    parts = line.split('\t') 
    chromosome = parts[0] 
    position = int(parts[2]) 
    stopcodon = parts[3] 
    flankseq = parts[4].strip() 
    flankseq = flankseq[3:6]+flankseq[9:12] 
    if flankseq in sequenceCount: 
     sequenceCount[flankseq] += 1  
for s in sequenceCount: 
    outfile.write(stopcodon+'\t'+s+'\t'+str(sequenceCount[s])+'\n') 

回答

2

您的outfile.write發生在for line in infile循環之外,因此stopcodon的值始終爲您在輸入文本文件的最後一行中的值。

如果您嘗試將序列計數與終止密碼子和側翼序列相關聯,則需要將這兩個變量用作關鍵字。如果您提前不知道所有終止密碼子,則無法使用「多重嵌套for循環」方法將sequenceCount的值初始化爲0,因此您應該使用defaultdict。

from collections import defaultdict 
sequenceCount = defaultdict(int) 

infile = open('flanking seqs.txt','rU') 
outfile = open('context resultsNEW.txt','w') 

for line in infile: 
    parts = line.split('\t') 
    chromosome = parts[0] 
    position = int(parts[2]) 
    stopcodon = parts[3] 
    flankseq = parts[4].strip() 
    flankseq = flankseq[3:6]+flankseq[9:12] 
    sequenceCount[flankseq, stopcodon] += 1 
for key, value in sequenceCount.iteritems(): 
    flankseq, stopcodon = key 
    outfile.write(stopcodon+'\t'+s+'\t'+str(sequenceCount[s])+'\n') 
0

當你出示你的輸出,您打印從每行的文件讀取的最後一個stopcodon,不管什麼stopcodon值在前面的循環使用。也許你的sequenceCount字典需要stopcodonflankseq的組合索引?