【问题标题】:no output for existing substring in a main string python/biopython主字符串 python/biopython 中的现有子字符串没有输出
【发布时间】:2021-03-23 15:33:14
【问题描述】:

如果主字符串的特定部分存在子字符串,我会构建一个代码来打印字符串。我有一个如下文件,我从seq11_rv 创建了 5 个字母子串(5mers)。

>seq11_fw
TCAGATGTGTATAAGAGACAGTTATTAGCCGGTTCCAGGTATGCAGTATGAGAA
>seq11_rv
GAGATTATGTGGGAAAGTTCATGGAATCGAGCGGAGATGTGTATAAGAGACAGTGCCGCGCTTCACTAGAAGTCATACTGC

然后我对这些 5mers 进行反向补码并将它们附加到一个列表中。然后我查看了seq11_fw,如果位置[42:51](seq11_fw 中的GCAGTATGA)有任何列表项,则应打印确认。

为了便于理解seq11_rv 的最后 5mer 是 ACTGC,它的反向补码变为 GCAGT,如果您检查 seq11_fw[42:51],则此 GCAGT 存在于该位置内,但我知道没有得到任何输出。

任何帮助将不胜感激。

这是我的代码:

from Bio import SeqIO
from Bio.Seq import Seq

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])


        cnt=0
        if len(record.seq) == 54 and any(str(items) in str(record.seq[42:51]) for items in lst):
            cnt +=1

        if cnt == 1:
            print(record.id)
            print(record.seq)

    print(lst)

【问题讨论】:

  • print(cnt) 怎么说?
  • print(any(str(items) in str(record.seq[42:51]) for items in lst)) 说什么?

标签: python performance output biopython


【解决方案1】:

这似乎可行,问题在于您设置算法的方式:

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 列表

第二个检查列表中长度 = 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() 

【讨论】:

    猜你喜欢
    • 2018-04-02
    • 2012-09-14
    • 1970-01-01
    • 2014-08-07
    • 2021-03-30
    • 2013-03-22
    • 1970-01-01
    • 1970-01-01
    • 2015-05-19
    相关资源
    最近更新 更多