用re.search()对基因组tRNA序列进行排序

Question

使用Biopython为tRNAs导入fasta格式的基因组序列，我正在编写这个搜索脚本；当我搜索'AT‘时，在这里下载的文件中：http://gtrnadb.ucsc.edu/Hsapi19/hg19-tRNAs.fa；返回626个结果中的622个。在我想要排除的4项搜索中，“AT”序列实际上存在于tRNA 54中(如下所示)。为什么re.search()没有将此结果包含在622中？我不确定我是误用了re.search()还是字符串太长了？‘as’可以按该顺序出现字符串中的任何内容，如下所述：https://docs.python.org/2/library/re.html

    #!/usr/bin/env python                                                                                                                                                                                                                
    import sys                                                                                                                                                                                                                           
    import os                                                                                                                                                                                                                            
    import time                                                                                                                                                                                                                          
    import re                                                                                                                                                                                                                            


    from Bio import SeqIO                                                                                                                                                                                                                

    filename = 'hg19-tRNAs270.fa'                                                                                                                                                                                                        
    search_seq = sys.argv[1]                                                                                                                                                                                                             

    count = 0                                                                                                                                                                                                                            
    fails = 0                                                                                                                                                                                                                            
    finds = 0                                                                                                                                                                                                                            

    output_handle = open("output_%s.fasta" % search_seq, "w")                                                                                                                                                                            

    for record in SeqIO.parse(filename, "fasta"):                                                                                                                                                                                        
        count += 1                                                                                                                                                                                                                       

        # print("Record " + record.id + ", length " + str(len(record.seq)))                                                                                                                                                              
        # print record                                                                                                                                                                                                                   
        # print 'ID:',record.id                                                                                                                                                                                                          
        # print 'SEQ:',record.seq                                                                                                                                                                                                        
        # print 'DESC:',record.description                                                                                                                                                                                               
        # print 'LENGTH:',str(len(record.seq))                                                                                                                                                                                           

        str_record = record.format("fasta")                                                                                                                                                                                              
        result = re.search(search_seq,str_record)                                                                                                                                                                                        


        if result != None:                                                                                                                                                                                                               
            finds +=1                                                                                                                                                                                                                    
            SeqIO.write(record, output_handle, "fasta")                                                                                                                                                                                  
        else:                                                                                                                                                                                                                            
            print 'ID:',record.id                                                                                                                                                                                                        
            print 'SEQ:',record.seq                                                                                                                                                                                                      
            print 'DESC:',record.description                                                                                                                                                                                             
            print 'LENGTH:',str(len(record.seq))                                                                                                                                                                                         


    output_handle.close()                                                                                                                                                                                                                
    print("There were " + str(count) + " records in file " + filename)                                                                                                                                                                   
    print("There were " + str(finds) + " matches in the file " + filename + " for sequence: " + (search_seq))

-模式:编译；默认目录："~/Dropbox/shared/“--

汇编于星期五7月11日15:35:06开始。

./fasta_search.py AT

ID: complimentry_to_oligo_probe

索赔序号：

DESC: complimentry_to_oligo_probe TCAGTTGGTAGAGCGGAGGA

长度:0

ID: Homo_sapiens_chr1.trna25-AsnGTT

SEQ：

DESC:人智人_chr1.trna25-AsnGTT (148000805-148000878) Asn (GTT) 74 bp Sc: 72.61

长度: 74

ID:人智人_chr19.trna12-GlnTTG

SEQ：

DESC:人Hr19.trna12-GlnTTG (9150428-9150356) Gln (TTG) 73 bp Sc: 23.15

长度: 73

ID:智人_chr6.trna54-ThrCGT

SEQ: GGCCCTGTAGCTCAGCGGTTGGAGCGCTGGTCTCGTAAACCTAGGGGTCGTGAGTTCAA*AT*CTCACCAGGGCCT

DESC:智人_chr6.trna54-ThrCGT (27586135-27586208) Thr (CGT) 74 bp Sc: 63.29

长度: 74

在hg19-tRNAs270.fa文件中有626条记录

文件hg19-tRNAs270.fa中有622个匹配序列: AT

汇编已于7月11日(星期五) 15:35:06完成

Answer 1

可能是断线。FASTA记录通常限制在每一行60个字符(或者web告诉我的)，所以如果格式调用将匹配的部分转换为A\nT，那么它将与您的regexp不匹配。使用str_record = str(record.seq)代替(假设您真的只想搜索序列)。

Answer 2

变化

result = re.search(search_seq,str_record)   

至

result = re.search(search_seq,str_record, re.MULTILINE)

我猜record.format()为您的str_record字符串添加了一个换行符。