|
import os
def get_oneseq(chr, start, end): # chr1, 23, 65
f = open(chr+'.fa', 'r')
head = f.readline() # >chr1\n
firstseqline = f.readline() # taaccctaaccctaaccctaaccctaaccctaaccctaaccctaacccta\n
offset = len(head) # 6
linelen = len(firstseqline) # 51
startpos = offset + start % (linelen-1) + (start/(linelen-1))*linelen-1
endpos = offset + end % (linelen-1) + (end/(linelen-1))*linelen-1
f.seek(startpos)
seq = f.read(endpos-startpos+1)
sequence = seq.replace(os.linesep,'')
return sequence
def parsebed(file,header = 0):
f3 = open(bedfile,'r')
f4 = open(file+'_seq.txt','w')
if header == 1:
newlines = f3.read().split(os.linesep)[1:-1]
else:
newlines = f3.read().split(os.linesep)[:-1]
seq = []
for i,c in enumerate(newlines):
chr = c.split()[0]
start = c.split()[1]
end = c.split()[2]
name = i+1
sequence = get_oneseq(chr, int(start), int(end))
print >>f4, '>'+str(name)
print >>f4, sequence
f3.close()
f4.close()
if __name__ == '__main__':
file = 'test.txt'
parsebed(file,1)
|
