0

给定一个 FASTA 文本文件 (Rosalind_gc.txt),我应该检查每个 DNA 记录并确定鸟嘌呤-胞嘧啶 (GC) 含量的百分比 (%)。

这方面的例子是:

样本数据集:

>Rosalind_6404
CCTGCGGAAGATCGGCACTAGAATAGCCAGAACCGTTTCTCTGAGGCTTCCGGCCTTCCC
TCCCACTAATAATTCTGAGG    
>Rosalind_5959
CCATCGGTAGCGCATCCTTAGTCCAATTAAGTCCCTATCCAGGCGCTCCGCCGAAGGTCT
ATATCCATTTGTCAGCAGACACGC
>Rosalind_0808
CCACCCTCGTGGTATGGCTAGGCATTCAGGAACCGGAGAACGCTTCAGACCAGCCCGGAC
TGGGAACCTGCGGGCAGTAGGTGGAAT

样本输出:

罗莎琳德_0808 60.919540

所以基本上遍历每个字符串,计算 G/C 出现的次数,然后将总数除以每个字符串的长度。我的问题是学习如何识别代码中的中断(即 >Rosalind_6404 )。我想要一个不使用 Biopython 和 biopython 方法的代码示例。

4

3 回答 3

2

您可以逐行读取文件并将序列数据累积到以“>”开头的下一行(再加上一次文件末尾)

def getCount(seq):
    return seq.count("G")+seq.count("C") 

with open("input.txt","r") as file:
    sequence = ""
    name     = ""
    for line in file:
        line = line.strip()
        if not line.startswith(">"):
            sequence += line
            continue
        if name != "":
            print(name, 100*getCount(sequence)/len(sequence))
        name     = line[1:]
        sequence = ""
    print(name, 100*getCount(sequence)/len(sequence))

# Rosalind_6404 53.75
# Rosalind_5959 53.57142857142857
# Rosalind_0808 60.91954022988506
于 2019-05-30T21:33:03.750 回答
2

由于您正在寻找 Biopython 解决方案,这里有一个非常简单的解决方案:

from Bio import SeqIO
from Bio.SeqUtils import GC

for r in SeqIO.parse('Rosalind_gc.fa', 'fasta'):
    print(r.id, GC(r.seq))

输出:

Rosalind_6404 53.75
Rosalind_5959 53.57142857142857
Rosalind_0808 60.91954022988506
于 2019-05-31T11:45:06.867 回答
1

您可能希望尽可能多地使用预编译的 C 模块来解决性能问题。有一个使用正则表达式的解决方案:

seq = 'CCTGCGGAAGATCGGCACTAGAATAGCCAGAACCGTTTCTCTGAGGCTTCCGGCCTTCCCTCCCACTAATAATTCTGAGG'

import re
perc = re.subn(r'[GC]', '', seq) / len(seq)

并处理“>”行:

seq = []
name = ''

for line in open('Rosalind_gc.txt'):
    if not line.startswith('>'):
        seq.append(line.strip())
    else:
        if seq:
            seq = ''.join(seq)
            perc = re.subn(r'[GC]', '', seq) / len(seq)
            print('{} has GC percent: {}'.format(name, perc * 100))
            seq = []
        name = line.strip()
于 2019-05-30T21:35:26.380 回答