问题描述
如果主字符串的特定部分存在子字符串,我会构建一个代码来打印字符串。我有一个如下文件,我从 seq11_rv
创建了 5 个字母子字符串 (5mers)。
>seq11_fw
TCAGATGTGTATAAGAGACAGTTATTAGCCGGTTCCAGGTATGCAGTATGAGAA
>seq11_rv
GAGATTATGTGGGAAAGTTCATGGAATCGAGCGGAGATGTGTATAAGAGACAGTGCCGCGCTTCACTAGAAGTCATACTGC
然后我对这些 5mer 进行反向补充并将它们附加到列表中。然后我查看了 seq11_fw
,如果位置 [42:51]
(seq11_fw 中的 GCAGTATGA)具有列表中的任何项目,则应打印确认。
为了更容易理解 seq11_rv
的最后 5mer 是 ACTGC
,它的反向补码变成 GCAGT
,如果你检查 seq11_fw[42:51]
这个 {{ 1}} 存在于该位置,但我没有得到任何输出。
任何帮助将不胜感激。
这是我的代码:
GCAGT
解决方法
这似乎可行,问题在于您设置算法的方式:
from Bio import SeqIO
from Bio.Seq import Seq
file ='test.faa'
with open(file,'r') as f:
lst = []
for record in SeqIO.parse(f,'fasta'):
if len(record.seq) == 81:
for i in range(len(record.seq)):
kmer = str(record.seq[i:i + 5])
if len(kmer) == 5:
C_kmer = Seq(kmer).complement()
lst.append(C_kmer[::-1])
with open(file,'r') as f:
cnt=0
for record in SeqIO.parse(f,'fasta'):
if len(record.seq) == 54 and any(str(items) in str(record.seq[41:52]) for items in lst):
cnt +=1
if cnt == 1:
print(record.id)
print(record.seq)
cnt = 0
在这个版本中,你对你的输入文件进行了两次迭代:
第一次创建 5mers 列表
第二个检查列表中 lenght = 54 的 record.seq 的 [41:52]。
您需要将计数器 cnt 重置为零,否则将打印所有序列
如果您计划拥有包含多个序列的测试文件,我会尝试使用具有相同顺序的一个具有正向序列的输入文件和一个具有反向序列的输入文件。然后我会检查每对并行循环两个文件的夫妇。不知道怎么做,但这里有很多很好的例子:
How to iterate through two lists in parallel?
我相信您需要使用 SeqIO.parse 将两条记录读入两个列表,然后像这样处理它们:
from Bio import SeqIO
from Bio.Seq import Seq
file1 ='test_fw.faa'
file2 ='test_rv.faa'
record1_lst = []
record2_lst = []
with open(file1,'r') as f1:
for record in SeqIO.parse(f1,'fasta'):
record1_lst.append(record)
with open(file2,'r') as f2:
for record in SeqIO.parse(f2,'fasta'):
record2_lst.append(record)
for record_fw,record_rv in zip(record1_lst,record2_lst):
print(record_fw.id,record_rv.id)
或者像这样工作,记住在解析器迭代器清空后关闭文件:
from Bio import SeqIO
from Bio.Seq import Seq
file1 ='test_fw.faa'
file2 ='test_rv.faa'
f1 = open(file1,'r')
f2 = open(file2,'r')
record1 = SeqIO.parse(f1,'fasta')
record2 = SeqIO.parse(f2,'fasta')
for record_fw,record_rv in zip(record1,record2):
print(record_fw.id,record_rv.id)
f1.close()
f2.close()