先把所有的fasta 序列合并为一个文件

cat *.fasta > result_all_Sequence.fasta

筛选大于90bp的

然后查找其中的GCAT

from glob import glob
import os

os.chdir("D:\\")

file_in = open("result_all_Sequences.fasta", ‘r‘) #定义文件file_in，为打开文件result_all_Sequences.fasta


fa_Con = file_in.read()           #.read()是把文件的全部内容读进来

file_in.close()                   #.close()是把文件关闭

every_fas =  fa_Con.split(">")    #.split(">")是指以>为分隔符把字符串分割为列表，分割后的列表里面不会包含>,即分割后>消失

## 写入文件

out_file = open("res.fasta", ‘w‘) #  w是可写

for i in every_fas:
    if i != "":
        start = i.index("\n")
        if len(i[start:]) >= 90:
            out_file.write(">" + i)
out_file.close()


# 读取刚刚筛选过的文件
file_in = open("res.fasta", ‘r‘) #定义文件file_in，为打开文件result_all_Sequences.fasta


fa_Con = file_in.read()          #.read()是把文件的全部内容读进来

file_in.close()                  #.close()是把文件关闭

every_fas =  fa_Con.split(">")   #.split(">")是指以>为分隔符把字符串分割为列表，分割后的列表里面不会包含>,即分割后>消失


out_file = open("ress.fasta", ‘w‘)

## 查找每条序列的GCAT
for i in every_fas:
    if i != "":
        start = i.index("\n")    #.index("\n")是指显示\n的位置或下标
    
        #print(i[start:])        # i[start:]是从\n开始一直到最后

        seq_con = i[start:]
      
        st_1 = seq_con.find("GCAT")
        
        if st_1 != -1 and st_1 != len(seq_con) - 4:
            out_file.write(">" + i[:start+1]+seq_con[st_1+4:])

   
out_file.close()

本文出自 “R和Python应用” 博客，请务必保留此出处http://matrix6ro.blog.51cto.com/1746429/1885915

python根据序列长度筛选序列，查找序列中指定碱基序列