2017-05-27 120 views
-1

我在Python中同时计算核苷酸和序列时遇到了问题。这是fasta文件,我需要对核苷酸和序列进行计数。请参考下面应该怎样是asnwer:计算核苷酸和序列

>chr12_9180206_+:chr12_118582391_+:a1;2 total_counts: 115 Seed: 4 K: 20 length: 79 
TTGGTTTCGTGGTTTTGCAAAGTATTGGCCTCCACCGCTATGTCTGGCTGGTTTACGA 
GCAGGACAGGCCGCTAAAGTG 
>chr12_9180206_+:chr12_118582391_+:a2;2 total_counts: 135 Seed: 4 K: 20 length: 80 
CTAACCCCCTACTTCCCAGACAGCTGCTCGTACAGTTTGGGCACATAGTCATCCCACTCG 
GCCTGGTAACACGTGCCAGC 
>chr1_8969882_-:chr1_568670_-:a1;113 total_counts: 7600 Seed: 225 K: 20 length: 86 
CACTCATGAGCTGTCCCCACATTAGGCTTAAAAACAGATGCAATTCCCGGACGTCTAAAC 
CAAACCACTTTCACCGCCACACGACC 
>chr1_8969882_-:chr1_568670_-:a2;69 total_counts: 6987 Seed: 197 K: 20 length: 120 
TGAACCTACGACTACACCGACTACGGCGGACTAATCTTCAACTCCTACATACTTCCCCCA 
TTATTCCTAGAACCAGGCGACCTGCGACTCCTTGACGTTGACAATCGAGTAGTACTCCCG 

而且我的代码:

f= open("data/assembledSeqs.fa", 'r') 
texto =f.read() 
f.close() 
TotalA=0 
TotalC=0 
TotalG=0 
TotalT=0 
cont=0 

for linea in texto.split('\n'): 
    if linea.startswith('>'): 
    print ("Secuencia: %d") % cont 
    cont+=1 
else: 
    TotalA=linea.count('A') 
    TotalC=linea.count('C') 
    TotalG=linea.count('G') 
    TotalT=linea.count('T') 
    print("Numero de A's: %d")%TotalA 
    print("Numero de C's: %d")%TotalC 
    print("Numero de G's: %d")%TotalG 
    print("Numero de T's: %d")%TotalT 

它的回报:

Secuencia: 0 
Numero de A's: 8 
Numero de C's: 12 
Numero de G's: 16 
Numero de T's: 22 

Numero de A's: 6 
Numero de C's: 5 
Numero de G's: 8 
Numero de T's: 2 

Secuencia: 1 
Numero de A's: 13 
Numero de C's: 23 
Numero de G's: 10 
Numero de T's: 14 

Numero de A's: 4 
Numero de C's: 7 
Numero de G's: 6 
Numero de T's: 3 

Secuencia: 2 
Numero de A's: 19 
Numero de C's: 18 
Numero de G's: 10 
Numero de T's: 13 

Numero de A's: 8 
Numero de C's: 13 
Numero de G's: 2 
Numero de T's: 3 

Secuencia: 3 
Numero de A's: 17 
Numero de C's: 23 
Numero de G's: 7 
Numero de T's: 13 

Numero de A's: 14 
Numero de C's: 18 
Numero de G's: 13 
Numero de T's: 15 

而且我想:

Secuencia 0: 
Número de A's: 14 
Número de C's: 17 
Número de G's: 24 
Número de T's: 24 

Secuencia 1: 
Número de A's: 17 
Número de C's: 30 
Número de G's: 16 
Número de T's: 17 

Secuencia 2: 
Número de A's: 27 
Número de C's: 31 
Número de G's: 12 
Número de T's: 16 

Secuencia 3: 
Número de A's: 31 
Número de C's: 41 
Número de G's: 20 
Número de T's: 28 
+0

你向我求助几次 - 所以请如果我的解决方案工作,然后不删除/破坏你的问题,请考虑[接受](https://meta.stackexchange.com/questions/5234/how-does-accepting-an-answer-work)第e为你工作的答案(upvoting/accepting相当于在StackOverflow上说“谢谢”)。 – MSeifert

回答

0

我就用一个字符串,而不是一个文件,但在理论上,你可以只打开文件并读取其中的内容:f = open('your_filename').read()

f = """>chr12_9180206_+:chr12_118582391_+:a1;2 total_counts: 115 Seed: 4 K: 20 length: 79 
TTGGTTTCGTGGTTTTGCAAAGTATTGGCCTCCACCGCTATGTCTGGCTGGTTTACGA 
GCAGGACAGGCCGCTAAAGTG 
>chr12_9180206_+:chr12_118582391_+:a2;2 total_counts: 135 Seed: 4 K: 20 length: 80 
CTAACCCCCTACTTCCCAGACAGCTGCTCGTACAGTTTGGGCACATAGTCATCCCACTCG 
GCCTGGTAACACGTGCCAGC 
>chr1_8969882_-:chr1_568670_-:a1;113 total_counts: 7600 Seed: 225 K: 20 length: 86 
CACTCATGAGCTGTCCCCACATTAGGCTTAAAAACAGATGCAATTCCCGGACGTCTAAAC 
CAAACCACTTTCACCGCCACACGACC 
>chr1_8969882_-:chr1_568670_-:a2;69 total_counts: 6987 Seed: 197 K: 20 length: 120 
TGAACCTACGACTACACCGACTACGGCGGACTAATCTTCAACTCCTACATACTTCCCCCA 
TTATTCCTAGAACCAGGCGACCTGCGACTCCTTGACGTTGACAATCGAGTAGTACTCCCG""" 

的第一步是由新行分割:

lines = f.split('\n') 

然后它,直到你发现与'>chr'开始的行只是一个收集计数的事情:

from collections import Counter 
cnts = [] 
last = Counter() 
for line in lines: 
    if line.startswith('>chr'): 
     if last: # only do that for the second, third, etc. 
      cnts.append(last) 
      last = Counter() 
    else: 
     last = last + Counter(line) 
cnts.append(last) # also append the last counts 

然后将收集的计数是:

>>> print(cnts) 
[Counter({'G': 24, 'T': 24, 'C': 17, 'A': 14}), 
Counter({'C': 30, 'A': 17, 'T': 17, 'G': 16}), 
Counter({'C': 31, 'A': 27, 'T': 16, 'G': 12}), 
Counter({'C': 41, 'A': 31, 'T': 28, 'G': 20})] 

其他的一切只是格式化的问题:

for seqno, cnt in enumerate(cnts): 
    print('Secuencia {}:'.format(seqno)) 
    for amino in ('A', 'C', 'G', 'T'): 
     print("Número de {}'s: {}".format(amino, cnt[amino])) 

它打印:

Secuencia 0: 
Número de A's: 14 
Número de C's: 17 
Número de G's: 24 
Número de T's: 24 
Secuencia 1: 
Número de A's: 17 
Número de C's: 30 
Número de G's: 16 
Número de T's: 17 
Secuencia 2: 
Número de A's: 27 
Número de C's: 31 
Número de G's: 12 
Número de T's: 16 
Secuencia 3: 
Número de A's: 31 
Número de C's: 41 
Número de G's: 20 
Número de T's: 28 
0

你可以使用str.count()然后解析该行分开。

+0

它不起作用,对不起。 – Mikel