我正在尝试反向补充 fasta DNA 序列
I am trying to reverse complement a fasta DNA sequence
我一直在尝试反向互补fasta DNA序列。这是我的代码:
fastafile=open('sequence (3).fasta','r')
entries=[]
reverse=""
sequence=['A','T','G','C','N']
for line in fastafile:
if not line.startswith('>'):
line = line.split()
entries.append(line)
print entries
for index in range(0,len(entries[::-1])):
if index !=sequence:
print "this is not a valid nucleotide"
break
else:
if index=='A':
reverse+='T'
elif index=='T':
reverse+='A'
elif index=='C':
reverse+='G'
elif index=='G':
reverse+ 'C'
elif index=='N':
reverse+='N'
print reverse
每次我得到输出时,这都不是一个有效的核苷酸,即使我打印的条目显示它有顺序的项目。这是我打印条目时的输出示例;
[['GCTCCCCTGAGGTTCGGCACCCACACTCCCTTCCCAGGAGCTCGCGATGCAAGAGCCACAGTCAGAGCTC'], ['AATATCGACCCCCCTCTGAGCCAGGAGACATTTTCAGAATTGTGGAACCTGCTTCCTGAAAACAATGTTC'], ['TGTCTTCGGAGCTGTGCCCAGCAGTGGATGAGCTGCTGCTCCCAGAGAGCGTCGTGAACTGGCTAGACGA']
如何解决这个问题?我只想补充一点,大约 2 个月前我才开始认真地使用 python 编程,所以我仍在学习和改进。
谢谢!
你的循环语句是:
for index in range(0,len(entries[::-1])):
这将遍历条目的长度,即 0, 1, 2, 3, ..., len(entries)
。
当你做 if index != sequence
时,你实际上是在比较一个整数和一个列表,比如 if 3 != ['A', 'C', 'T', 'G']
。我假设你可以看到这没有意义。您可能想要做的是查看序列中的核苷酸是否是有效核苷酸,因此它在 sequence
列表中。你可以这样做:
if entries[::-1][index] in sequence # Will be true if the nucleotide at entries[::-1][index] is inside sequence
说两句:
第一个,你不必范围到len(entries[::-1])
,它和len(entries)
一样
其次,更重要的是,有一个专门为生物信息学构建的实际模块。它被称为Biopython。它具有特殊的对象和功能。例如,您的问题可以按如下方式解决:
-
from Bio.Seq import Seq
dna = Seq("ATGGCCATTGTAATGGGCCGCTGAAAGGGTGCCCGATAG")
print dna.reverse_complement()
输出: CTATCGGGCACCCTTTCAGCGGCCCATTACAATGGCCAT
我一直在尝试反向互补fasta DNA序列。这是我的代码:
fastafile=open('sequence (3).fasta','r')
entries=[]
reverse=""
sequence=['A','T','G','C','N']
for line in fastafile:
if not line.startswith('>'):
line = line.split()
entries.append(line)
print entries
for index in range(0,len(entries[::-1])):
if index !=sequence:
print "this is not a valid nucleotide"
break
else:
if index=='A':
reverse+='T'
elif index=='T':
reverse+='A'
elif index=='C':
reverse+='G'
elif index=='G':
reverse+ 'C'
elif index=='N':
reverse+='N'
print reverse
每次我得到输出时,这都不是一个有效的核苷酸,即使我打印的条目显示它有顺序的项目。这是我打印条目时的输出示例;
[['GCTCCCCTGAGGTTCGGCACCCACACTCCCTTCCCAGGAGCTCGCGATGCAAGAGCCACAGTCAGAGCTC'], ['AATATCGACCCCCCTCTGAGCCAGGAGACATTTTCAGAATTGTGGAACCTGCTTCCTGAAAACAATGTTC'], ['TGTCTTCGGAGCTGTGCCCAGCAGTGGATGAGCTGCTGCTCCCAGAGAGCGTCGTGAACTGGCTAGACGA']
如何解决这个问题?我只想补充一点,大约 2 个月前我才开始认真地使用 python 编程,所以我仍在学习和改进。 谢谢!
你的循环语句是:
for index in range(0,len(entries[::-1])):
这将遍历条目的长度,即 0, 1, 2, 3, ..., len(entries)
。
当你做 if index != sequence
时,你实际上是在比较一个整数和一个列表,比如 if 3 != ['A', 'C', 'T', 'G']
。我假设你可以看到这没有意义。您可能想要做的是查看序列中的核苷酸是否是有效核苷酸,因此它在 sequence
列表中。你可以这样做:
if entries[::-1][index] in sequence # Will be true if the nucleotide at entries[::-1][index] is inside sequence
说两句:
第一个,你不必范围到
len(entries[::-1])
,它和len(entries)
一样
其次,更重要的是,有一个专门为生物信息学构建的实际模块。它被称为Biopython。它具有特殊的对象和功能。例如,您的问题可以按如下方式解决:
-
from Bio.Seq import Seq
dna = Seq("ATGGCCATTGTAATGGGCCGCTGAAAGGGTGCCCGATAG")
print dna.reverse_complement()
输出: CTATCGGGCACCCTTTCAGCGGCCCATTACAATGGCCAT