annotate bedutil.py @ 9:4d3ac66875d2 draft default tip

Uploaded
author jjohnson
date Thu, 14 Dec 2017 13:35:00 -0500
parents 5c92d0be6514
children
Ignore whitespace changes - Everywhere: Within whitespace: At end of lines:
rev   line source
8
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
1 #!/usr/bin/env python
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
2 """
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
3 #
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
4 #------------------------------------------------------------------------------
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
5 # University of Minnesota
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
6 # Copyright 2016, Regents of the University of Minnesota
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
7 #------------------------------------------------------------------------------
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
8 # Author:
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
9 #
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
10 # James E Johnson
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
11 #
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
12 #------------------------------------------------------------------------------
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
13 """
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
14
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
15 import sys
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
16 from Bio.Seq import reverse_complement, translate
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
17
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
18
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
19 def bed_from_line(line,ensembl=False):
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
20 fields = line.rstrip('\r\n').split('\t')
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
21 if len(fields) < 12:
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
22 return None
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
23 (chrom, chromStart, chromEnd, name, score, strand,
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
24 thickStart, thickEnd, itemRgb,
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
25 blockCount, blockSizes, blockStarts) = fields[0:12]
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
26 bed_entry = BedEntry(chrom=chrom, chromStart=chromStart, chromEnd=chromEnd,
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
27 name=name, score=score, strand=strand,
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
28 thickStart=thickStart, thickEnd=thickEnd,
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
29 itemRgb=itemRgb,
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
30 blockCount=blockCount,
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
31 blockSizes=blockSizes.rstrip(','),
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
32 blockStarts=blockStarts.rstrip(','))
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
33 if ensembl and len(fields) >= 20:
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
34 bed_entry.second_name = fields[12]
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
35 bed_entry.cds_start_status = fields[13]
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
36 bed_entry.cds_end_status = fields[14]
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
37 bed_entry.exon_frames = fields[15].rstrip(',')
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
38 bed_entry.biotype = fields[16]
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
39 bed_entry.gene_name = fields[17]
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
40 bed_entry.second_gene_name = fields[18]
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
41 bed_entry.gene_type = fields[19]
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
42 return bed_entry
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
43
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
44
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
45
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
46 class BedEntry( object ):
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
47 def __init__(self, chrom=None, chromStart=None, chromEnd=None,
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
48 name=None, score=None, strand=None,
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
49 thickStart=None, thickEnd=None, itemRgb=None,
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
50 blockCount=None, blockSizes=None, blockStarts=None):
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
51 self.chrom = chrom
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
52 self.chromStart = int(chromStart)
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
53 self.chromEnd = int(chromEnd)
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
54 self.name = name
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
55 self.score = int(score) if score is not None else 0
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
56 self.strand = '-' if str(strand).startswith('-') else '+'
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
57 self.thickStart = int(thickStart) if thickStart else self.chromStart
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
58 self.thickEnd = int(thickEnd) if thickEnd else self.chromEnd
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
59 self.itemRgb = str(itemRgb) if itemRgb is not None else r'100,100,100'
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
60 self.blockCount = int(blockCount)
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
61 if isinstance(blockSizes, str) or isinstance(blockSizes, unicode):
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
62 self.blockSizes = [int(x) for x in blockSizes.split(',')]
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
63 elif isinstance(blockSizes, list):
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
64 self.blockSizes = [int(x) for x in blockSizes]
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
65 else:
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
66 self.blockSizes = blockSizes
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
67 if isinstance(blockStarts, str) or isinstance(blockSizes, unicode):
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
68 self.blockStarts = [int(x) for x in blockStarts.split(',')]
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
69 elif isinstance(blockStarts, list):
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
70 self.blockStarts = [int(x) for x in blockStarts]
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
71 else:
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
72 self.blockStarts = blockStarts
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
73 self.second_name = None
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
74 self.cds_start_status = None
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
75 self.cds_end_status = None
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
76 self.exon_frames = None
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
77 self.biotype = None
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
78 self.gene_name = None
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
79 self.second_gene_name = None
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
80 self.gene_type = None
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
81 self.seq = None
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
82 self.cdna = None
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
83 self.pep = None
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
84 # T26C
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
85 self.aa_change = []
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
86 # p.Trp26Cys g.<pos><ref>><alt> # g.1304573A>G
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
87 self.variants = []
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
88
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
89
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
90 def __str__(self):
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
91 return '%s\t%d\t%d\t%s\t%d\t%s\t%d\t%d\t%s\t%d\t%s\t%s' % (
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
92 self.chrom, self.chromStart, self.chromEnd,
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
93 self.name, self.score, self.strand,
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
94 self.thickStart, self.thickEnd, str(self.itemRgb), self.blockCount,
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
95 ','.join([str(x) for x in self.blockSizes]),
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
96 ','.join([str(x) for x in self.blockStarts]))
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
97
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
98
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
99 def get_splice_junctions(self):
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
100 splice_juncs = []
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
101 for i in range(self.blockCount - 1):
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
102 splice_junc = "%s:%d_%d" % (self.chrom, self.chromStart + self.blockSizes[i], self.chromStart + self.blockStarts[i+1])
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
103 splice_juncs.append(splice_junc)
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
104 return splice_juncs
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
105
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
106
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
107 def get_exon_seqs(self):
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
108 if not self.seq:
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
109 return None
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
110 exons = []
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
111 for i in range(self.blockCount):
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
112 # splice_junc = "%s:%d_%d" % (self.chrom, self.chromStart + self.blockSizes[i], self.chromStart + self.blockStarts[i+1])
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
113 exons.append(self.seq[self.blockStarts[i]:self.blockStarts[i] + self.blockSizes[i]])
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
114 if self.strand == '-': #reverse complement
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
115 exons.reverse()
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
116 for i,s in enumerate(exons):
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
117 exons[i] = reverse_complement(s)
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
118 return exons
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
119
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
120
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
121 def get_spliced_seq(self,strand=None):
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
122 if not self.seq:
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
123 return None
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
124 seq = ''.join(self.get_exon_seqs())
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
125 if strand and self.strand != strand:
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
126 seq = reverse_complement(seq)
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
127 return seq
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
128
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
129
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
130 def get_cdna(self):
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
131 if not self.cdna:
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
132 self.cdna = self.get_spliced_seq()
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
133 return self.cdna
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
134
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
135
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
136 def get_cds(self):
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
137 cdna = self.get_cdna()
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
138 if cdna:
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
139 if self.chromStart == self.thickStart and self.chromEnd == self.thickEnd:
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
140 return cdna
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
141 pos = [self.offset_of_pos(self.thickStart),self.offset_of_pos(self.thickEnd)]
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
142 if 0 <= min(pos) <= max(pos) <= len(cdna):
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
143 return cdna[min(pos):max(pos)]
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
144 return None
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
145
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
146
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
147 def get_cigar(self):
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
148 cigar = ''
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
149 md = ''
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
150 r = range(self.blockCount)
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
151 rev = self.strand == '-'
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
152 ## if rev: r.reverse()
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
153 xl = None
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
154 for x in r:
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
155 if xl is not None:
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
156 intronSize = abs(self.blockStarts[x] - self.blockSizes[xl] - self.blockStarts[xl])
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
157 cigar += '%dN' % intronSize
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
158 cigar += '%dM' % self.blockSizes[x]
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
159 xl = x
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
160 return cigar
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
161
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
162
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
163 def get_cigar_md(self):
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
164 cigar = ''
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
165 md = ''
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
166 r = range(self.blockCount)
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
167 rev = self.strand == '-'
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
168 ## if rev: r.reverse()
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
169 xl = None
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
170 for x in r:
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
171 if xl is not None:
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
172 intronSize = abs(self.blockStarts[x] - self.blockSizes[xl] - self.blockStarts[xl])
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
173 cigar += '%dN' % intronSize
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
174 cigar += '%dM' % self.blockSizes[x]
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
175 xl = x
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
176 md = '%d' % sum(self.blockSizes)
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
177 return (cigar,md)
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
178
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
179
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
180 def get_translation(self,sequence=None):
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
181 translation = None
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
182 seq = sequence if sequence else self.get_spliced_seq()
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
183 if seq:
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
184 seqlen = len(seq) / 3 * 3;
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
185 if seqlen >= 3:
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
186 translation = translate(seq[:seqlen])
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
187 return translation
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
188
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
189
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
190 def get_translations(self):
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
191 translations = []
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
192 seq = self.get_spliced_seq()
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
193 if seq:
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
194 for i in range(3):
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
195 translation = self.get_translation(sequence=seq[i:])
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
196 if translation:
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
197 translations.append(translation)
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
198 return translations
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
199
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
200 def pos_of_na(self, offset):
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
201 if offset is not None and 0 <= offset < sum(self.blockSizes):
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
202 r = range(self.blockCount)
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
203 rev = self.strand == '-'
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
204 if rev:
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
205 r.reverse()
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
206 nlen = 0
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
207 for x in r:
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
208 if offset < nlen + self.blockSizes[x]:
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
209 if rev:
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
210 return self.chromStart + self.blockStarts[x] + self.blockSizes[x] - (offset - nlen)
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
211 else:
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
212 return self.chromStart + self.blockStarts[x] + (offset - nlen)
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
213 nlen += self.blockSizes[x]
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
214 return None
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
215
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
216 def offset_of_pos(self, pos):
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
217 if not self.chromStart <= pos < self.chromEnd:
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
218 return -1
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
219 r = range(self.blockCount)
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
220 rev = self.strand == '-'
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
221 if rev:
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
222 r.reverse()
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
223 nlen = 0
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
224 for x in r:
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
225 bStart = self.chromStart + self.blockStarts[x]
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
226 bEnd = bStart + self.blockSizes[x]
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
227 ## print >> sys.stdout, "offset_of_pos %d %d %d %d" % (bStart,pos,bEnd,nlen)
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
228 if bStart <= pos < bEnd:
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
229 return nlen + (bEnd - pos if rev else pos - bStart)
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
230 nlen += self.blockSizes[x]
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
231
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
232 def apply_variant(self,pos,ref,alt):
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
233 pos = int(pos)
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
234 if not ref or not alt:
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
235 print >> sys.stderr, "variant requires ref and alt sequences"
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
236 return
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
237 if not self.chromStart <= pos <= self.chromEnd:
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
238 print >> sys.stderr, "variant not in entry %s: %s %d < %d < %d" % (self.name,self.strand, self.chromStart,pos,self.chromEnd)
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
239 print >> sys.stderr, "%s" % str(self)
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
240 return
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
241 if len(ref) != len(alt):
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
242 print >> sys.stderr, "variant only works for snp: %s %s" % (ref,alt)
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
243 return
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
244 if not self.seq:
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
245 print >> sys.stderr, "variant entry %s has no seq" % self.name
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
246 return
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
247 """
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
248 if self.strand == '-':
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
249 ref = reverse_complement(ref)
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
250 alt = reverse_complement(alt)
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
251 """
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
252 bases = list(self.seq)
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
253 offset = pos - self.chromStart
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
254 for i in range(len(ref)):
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
255 # offset = self.offset_of_pos(pos+i)
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
256 if offset is not None:
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
257 bases[offset+i] = alt[i]
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
258 else:
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
259 print >> sys.stderr, "variant offset %s: %s %d < %d < %d" % (self.name,self.strand,self.chromStart,pos+1,self.chromEnd)
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
260 print >> sys.stderr, "%s" % str(self)
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
261 self.seq = ''.join(bases)
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
262 self.variants.append("g.%d%s>%s" % (pos+1,ref,alt))
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
263
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
264 def get_variant_bed(self,pos,ref,alt):
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
265 pos = int(pos)
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
266 if not ref or not alt:
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
267 print >> sys.stderr, "variant requires ref and alt sequences"
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
268 return None
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
269 if not self.chromStart <= pos <= self.chromEnd:
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
270 print >> sys.stderr, "variant not in entry %s: %s %d < %d < %d" % (self.name,self.strand, self.chromStart,pos,self.chromEnd)
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
271 print >> sys.stderr, "%s" % str(self)
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
272 return None
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
273 if not self.seq:
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
274 print >> sys.stderr, "variant entry %s has no seq" % self.name
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
275 return None
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
276 tbed = BedEntry(chrom = self.chrom, chromStart = self.chromStart, chromEnd = self.chromEnd,
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
277 name = self.name, score = self.score, strand = self.strand,
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
278 thickStart = self.chromStart, thickEnd = self.chromEnd, itemRgb = self.itemRgb,
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
279 blockCount = self.blockCount, blockSizes = self.blockSizes, blockStarts = self.blockStarts)
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
280 bases = list(self.seq)
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
281 offset = pos - self.chromStart
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
282 tbed.seq = ''.join(bases[:offset] + list(alt) + bases[offset+len(ref):])
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
283 if len(ref) != len(alt):
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
284 diff = len(alt) - len(ref)
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
285 rEnd = pos + len(ref)
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
286 aEnd = pos + len(alt)
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
287 #need to adjust blocks
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
288 # change spans blocks,
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
289 for x in range(tbed.blockCount):
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
290 bStart = tbed.chromStart + tbed.blockStarts[x]
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
291 bEnd = bStart + tbed.blockSizes[x]
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
292 # change within a block or extends (last block), adjust blocksize
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
293 # seq: GGGcatGGG
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
294 # ref c alt tag: GGGtagatGGG
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
295 # ref cat alt a: GGGaGGG
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
296 if bStart <= pos < rEnd < bEnd:
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
297 tbed.blockSizes[x] += diff
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
298 return tbed
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
299
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
300 ## (start,end)
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
301 def get_subrange(self,tstart,tstop,debug=False):
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
302 chromStart = self.chromStart
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
303 chromEnd = self.chromEnd
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
304 if debug:
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
305 print >> sys.stderr, "%s" % (str(self))
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
306 r = range(self.blockCount)
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
307 if self.strand == '-':
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
308 r.reverse()
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
309 bStart = 0
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
310 bEnd = 0
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
311 for x in r:
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
312 bEnd = bStart + self.blockSizes[x]
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
313 if bStart <= tstart < bEnd:
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
314 if self.strand == '+':
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
315 chromStart = self.chromStart + self.blockStarts[x] +\
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
316 (tstart - bStart)
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
317 else:
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
318 chromEnd = self.chromStart + self.blockStarts[x] +\
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
319 self.blockSizes[x] - (tstart - bStart)
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
320 if bStart <= tstop < bEnd:
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
321 if self.strand == '+':
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
322 chromEnd = self.chromStart + self.blockStarts[x] +\
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
323 (tstop - bStart)
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
324 else:
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
325 chromStart = self.chromStart + self.blockStarts[x] +\
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
326 self.blockSizes[x] - (tstop - bStart)
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
327 if debug:
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
328 print >> sys.stderr,\
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
329 "%3d %s\t%d\t%d\t%d\t%d\t%d\t%d"\
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
330 % (x, self.strand, bStart, bEnd,
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
331 tstart, tstop, chromStart, chromEnd)
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
332 bStart += self.blockSizes[x]
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
333 return(chromStart, chromEnd)
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
334
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
335
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
336 # get the blocks for sub range
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
337 def get_blocks(self, chromStart, chromEnd):
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
338 tblockCount = 0
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
339 tblockSizes = []
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
340 tblockStarts = []
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
341 for x in range(self.blockCount):
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
342 bStart = self.chromStart + self.blockStarts[x]
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
343 bEnd = bStart + self.blockSizes[x]
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
344 if bStart > chromEnd:
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
345 break
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
346 if bEnd < chromStart:
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
347 continue
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
348 cStart = max(chromStart, bStart)
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
349 tblockStarts.append(cStart - chromStart)
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
350 tblockSizes.append(min(chromEnd, bEnd) - cStart)
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
351 tblockCount += 1
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
352 return (tblockCount, tblockSizes, tblockStarts)
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
353
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
354
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
355 def trim(self, tstart, tstop, debug=False):
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
356 (tchromStart, tchromEnd) =\
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
357 self.get_subrange(tstart, tstop, debug=debug)
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
358 (tblockCount, tblockSizes, tblockStarts) =\
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
359 self.get_blocks(tchromStart, tchromEnd)
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
360 tbed = BedEntry(
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
361 chrom=self.chrom, chromStart=tchromStart, chromEnd=tchromEnd,
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
362 name=self.name, score=self.score, strand=self.strand,
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
363 thickStart=tchromStart, thickEnd=tchromEnd, itemRgb=self.itemRgb,
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
364 blockCount=tblockCount,
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
365 blockSizes=tblockSizes, blockStarts=tblockStarts)
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
366 if self.seq:
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
367 ts = tchromStart-self.chromStart
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
368 te = tchromEnd - tchromStart + ts
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
369 tbed.seq = self.seq[ts:te]
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
370 return tbed
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
371
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
372
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
373 def get_filtered_translations(self,untrimmed=False,filtering=True,ignore_left_bp=0,ignore_right_bp=0,debug=False):
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
374 translations = [None,None,None]
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
375 seq = self.get_spliced_seq()
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
376 ignore = (ignore_left_bp if self.strand == '+' else ignore_right_bp) / 3
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
377 block_sum = sum(self.blockSizes)
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
378 exon_sizes = [x for x in self.blockSizes]
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
379 if self.strand == '-':
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
380 exon_sizes.reverse()
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
381 splice_sites = [sum(exon_sizes[:x]) / 3 for x in range(1,len(exon_sizes))]
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
382 if debug:
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
383 print >> sys.stderr, "splice_sites: %s" % splice_sites
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
384 junc = splice_sites[0] if len(splice_sites) > 0 else exon_sizes[0]
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
385 if seq:
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
386 for i in range(3):
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
387 translation = self.get_translation(sequence=seq[i:])
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
388 if translation:
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
389 tstart = 0
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
390 tstop = len(translation)
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
391 offset = (block_sum - i) % 3
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
392 if debug:
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
393 print >> sys.stderr, "frame: %d\ttstart: %d tstop: %d offset: %d\t%s" % (i,tstart,tstop,offset,translation)
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
394 if not untrimmed:
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
395 tstart = translation.rfind('*',0,junc) + 1
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
396 stop = translation.find('*',junc)
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
397 tstop = stop if stop >= 0 else len(translation)
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
398 offset = (block_sum - i) % 3
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
399 trimmed = translation[tstart:tstop]
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
400 if debug:
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
401 print >> sys.stderr, "frame: %d\ttstart: %d tstop: %d offset: %d\t%s" % (i,tstart,tstop,offset,trimmed)
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
402 if filtering and tstart > ignore:
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
403 continue
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
404 #get genomic locations for start and end
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
405 if self.strand == '+':
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
406 chromStart = self.chromStart + i + (tstart * 3)
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
407 chromEnd = self.chromEnd - offset - (len(translation) - tstop) * 3
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
408 else:
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
409 chromStart = self.chromStart + offset + (len(translation) - tstop) * 3
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
410 chromEnd = self.chromEnd - i - (tstart * 3)
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
411 #get the blocks for this translation
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
412 (tblockCount,tblockSizes,tblockStarts) = self.get_blocks(chromStart,chromEnd)
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
413 translations[i] = (chromStart,chromEnd,trimmed,tblockCount,tblockSizes,tblockStarts)
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
414 if debug:
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
415 print >> sys.stderr, "tblockCount: %d tblockStarts: %s tblockSizes: %s" % (tblockCount,tblockStarts,tblockSizes)
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
416 # translations[i] = (chromStart,chromEnd,trimmed,tblockCount,tblockSizes,tblockStarts)
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
417 return translations
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
418
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
419
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
420 def get_seq_id(self,seqtype='unk:unk',reference='',frame=None):
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
421 ## Ensembl fasta ID format
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
422 # >ID SEQTYPE:STATUS LOCATION GENE TRANSCRIPT
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
423 # >ENSP00000328693 pep:splice chromosome:NCBI35:1:904515:910768:1 gene:ENSG00000158815:transcript:ENST00000328693 gene_biotype:protein_coding transcript_biotype:protein_coding
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
424 frame_name = ''
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
425 chromStart = self.chromStart
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
426 chromEnd = self.chromEnd
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
427 strand = 1 if self.strand == '+' else -1
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
428 if frame != None:
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
429 block_sum = sum(self.blockSizes)
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
430 offset = (block_sum - frame) % 3
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
431 frame_name = '_' + str(frame + 1)
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
432 if self.strand == '+':
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
433 chromStart += frame
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
434 chromEnd -= offset
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
435 else:
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
436 chromStart += offset
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
437 chromEnd -= frame
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
438 location = "chromosome:%s:%s:%s:%s:%s" % (reference,self.chrom,chromStart,chromEnd,strand)
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
439 seq_id = "%s%s %s %s" % (self.name,frame_name,seqtype,location)
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
440 return seq_id
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
441
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
442
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
443 def get_line(self, start_offset = 0, end_offset = 0):
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
444 if start_offset or end_offset:
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
445 s_offset = start_offset if start_offset else 0
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
446 e_offset = end_offset if end_offset else 0
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
447 if s_offset > self.chromStart:
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
448 s_offset = self.chromStart
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
449 chrStart = self.chromStart - s_offset
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
450 chrEnd = self.chromEnd + e_offset
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
451 blkSizes = self.blockSizes
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
452 blkSizes[0] += s_offset
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
453 blkSizes[-1] += e_offset
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
454 blkStarts = self.blockStarts
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
455 for i in range(1,self.blockCount):
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
456 blkStarts[i] += s_offset
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
457 items = [str(x) for x in [self.chrom,chrStart,chrEnd,self.name,self.score,self.strand,self.thickStart,self.thickEnd,self.itemRgb,self.blockCount,','.join([str(x) for x in blkSizes]),','.join([str(x) for x in blkStarts])]]
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
458 return '\t'.join(items) + '\n'
5c92d0be6514 Uploaded
jjohnson
parents:
diff changeset
459 return self.line