Problem
The GC-content of a DNA string is given by the percentage of symbols in the string that are ‘C’ or ‘G’. For example, the GC-content of “AGCTATAG” is 37.5%. Note that the reverse complement of any DNA string has the same GC-content.
DNA strings must be labeled when they are consolidated into a database. A commonly used method of string labeling is called FASTA format. In this format, the string is introduced by a line that begins with ‘>’, followed by some labeling information. Subsequent lines contain the string itself; the first line to begin with ‘>’ indicates the label of the next string.
In Rosalind’s implementation, a string in FASTA format will be labeled by the ID “Rosalind_xxxx”, where “xxxx” denotes a four-digit code between 0000 and 9999.
Given: At most 10 DNA strings in FASTA format (of length at most 1 kbp each).
Return: The ID of the string having the highest GC-content, followed by the GC-content of that string. Rosalind allows for a default error of 0.001 in all decimal answers unless otherwise stated; please see the note on absolute error below.
Sample Dataset
>Rosalind_6404
CCTGCGGAAGATCGGCACTAGAATAGCCAGAACCGTTTCTCTGAGGCTTCCGGCCTTCCC
TCCCACTAATAATTCTGAGG
>Rosalind_5959
CCATCGGTAGCGCATCCTTAGTCCAATTAAGTCCCTATCCAGGCGCTCCGCCGAAGGTCT
ATATCCATTTGTCAGCAGACACGC
>Rosalind_0808
CCACCCTCGTGGTATGGCTAGGCATTCAGGAACCGGAGAACGCTTCAGACCAGCCCGGAC
TGGGAACCTGCGGGCAGTAGGTGGAAT
Sample Output
Rosalind_0808
60.919540
Solution
#!/usr/bin/env python3
f=open("DATA/rosalind_gc.txt", "r")
seq = f.readlines()
gc = 0
n = 0
max_gc = 0
for s in seq:
s = s.strip()
if s[0] == '>':
if n == 0:
desc = s[1:]
continue
else:
if gc*100/n > max_gc:
max_gc = gc*100/n
max_desc = desc
n = 0
gc = 0
desc = s[1:]
else:
n += len(s)
gc += s.count('G')
gc += s.count('C')
if max_gc < gc*100/n:
max_gc = gc*100/n
max_desc = desc
print(max_desc)
print(max_gc)
这道题第一次出现了fasta文件，直接解析文件算gc含量，当我们第二次看到fasta文件的时候，就该写函数来解析了，而不是像上面这样子。预先上biopython的SeqIO来解析fasta，下次我们将自己写一个解析fasta的函数。
#!/usr/bin/env python3
from Bio import SeqIO
fasta_sequences = SeqIO.parse("DATA/rosalind_gc.txt",'fasta')
def gc(seq):
return (seq.count('G') + seq.count('C'))/len(seq)
id_vector = []
gc_vector = []
for fasta in fasta_sequences:
id_vector.append(fasta.id)
gc_vector.append(gc(str(fasta.seq)))
max_gc = max(gc_vector)
max_gc_idx = gc_vector.index(max_gc)
max_gc_id = id_vector[max_gc_idx]
print(max_gc_id)
print(max_gc * 100)