Back to index

python-biopython  1.60
Public Member Functions | Public Attributes
nextorf.NextOrf Class Reference

List of all members.

Public Member Functions

def __init__
def ReadFile
def ToFasta
def Gc
def Gc2
def GetOrfCoordinates
def GetCDS
def Output

Public Attributes

 options
 file
 genetic_code
 table
 counter
 header
 seq
 length
 rseq

Detailed Description

Definition at line 47 of file nextorf.py.


Constructor & Destructor Documentation

def nextorf.NextOrf.__init__ (   self,
  file,
  options 
)

Definition at line 48 of file nextorf.py.

00048 
00049     def __init__(self, file, options):
00050         self.options = options
00051         self.file = file
00052         self.genetic_code = int(self.options['table'])
00053         self.table = makeTableX(CodonTable.ambiguous_dna_by_id[self.genetic_code])
00054         self.counter = 0
00055         self.ReadFile()
        

Member Function Documentation

def nextorf.NextOrf.Gc (   self,
  seq 
)

Definition at line 85 of file nextorf.py.

00085 
00086     def Gc(self, seq):
00087        d = {}
00088        for nt in 'ATGC':
00089           d[nt] = seq.count(nt)
00090        gc = d['G'] + d['C']
00091        if gc == 0: return 0
00092        return round(gc*100.0/(d['A'] +d['T'] + gc),1)

def nextorf.NextOrf.Gc2 (   self,
  seq 
)

Definition at line 93 of file nextorf.py.

00093 
00094     def Gc2(self,seq):
00095        l = len(seq)
00096        d= {}
00097        for nt in ['A','T','G','C']:
00098           d[nt] = [0,0,0]
00099           
00100        for i in range(0,l,3):
00101           codon = seq[i:i+3]
00102           if len(codon) <3: codon = codon + '  '
00103           for pos in range(0,3):
00104              for nt in ['A','T','G','C']:
00105                 if codon[pos] == nt: d[nt][pos] = d[nt][pos] +1
00106 
00107        gc = {}
00108        gcall = 0
00109        nall = 0
00110        for i in range(0,3):
00111           try:
00112              n = d['G'][i] + d['C'][i] +d['T'][i] + d['A'][i]
00113              gc[i] = (d['G'][i] + d['C'][i])*100.0/n
00114           except:
00115              gc[i] = 0
00116 
00117           gcall = gcall + d['G'][i] + d['C'][i]
00118           nall = nall + n
00119 
00120        gcall = 100.0*gcall/nall
00121        res = '%.1f%%, %.1f%%, %.1f%%, %.1f%%' % (gcall, gc[0], gc[1], gc[2])
00122        return res

Here is the caller graph for this function:

def nextorf.NextOrf.GetCDS (   self,
  seq,
  strand = 1 
)

Definition at line 143 of file nextorf.py.

00143 
00144     def GetCDS(self, seq, strand = 1):
00145         frame_coordinates = self.GetOrfCoordinates(seq)
00146         START, STOP = 1,0
00147         so = self.options
00148         nostart = so['nostart']
00149         minlength, maxlength = int(so['minlength']), int(so['maxlength'])
00150         CDS = []
00151         f = 0
00152         for frame in frame_coordinates:
00153             f+=1
00154             start_site = 0
00155             if nostart == '1': start_site = 1
00156             frame.append((self.length, 0, 'XXX'))
00157             for pos, codon_type, codon in frame:
00158                 if codon_type == START:
00159                     if start_site == 0: start_site = pos
00160                 elif codon_type == STOP:
00161                     if start_site == 0: continue
00162 #                    if codon == 'XXX': print 'do something'
00163                     stop = pos + 2
00164 #                    print stop
00165                     length = stop - start_site +1
00166                     if length >= minlength and length <= maxlength:
00167                         if nostart == '1' and start_site == 1:
00168                            start_site = start_site + f - 1
00169                         if codon == 'XXX': stop = start_site + 3*((int((stop-1)-start_site)/3))
00170                         s = seq[start_site -1 : stop]
00171                         CDS.append((start_site, stop, length, s, strand*f))
00172                         start_site = 0
00173                         if nostart == '1': start_site = stop + 1
00174                     elif length < minlength or length > maxlength:
00175                        start_site = 0
00176                        if nostart == '1': start_site = stop + 1
00177                     del stop   
00178         return CDS

Here is the call graph for this function:

def nextorf.NextOrf.GetOrfCoordinates (   self,
  seq 
)

Definition at line 123 of file nextorf.py.

00123 
00124     def GetOrfCoordinates(self, seq):
00125         s = seq.data
00126         letters = []
00127         table = self.table
00128         get = self.table.forward_table.get
00129         n = len(seq)
00130         start_codons = self.table.start_codons
00131         stop_codons = self.table.stop_codons
00132 #        print 'Start codons', start_codons
00133 #        print 'Stop codons', stop_codons
00134         frame_coordinates = []
00135         for frame in range(0,3):
00136             coordinates = []
00137             for i in range(0+frame, n-n%3, 3):
00138                 codon = s[i:i+3]
00139                 if codon in start_codons: coordinates.append((i+1,1,codon))
00140                 elif codon in stop_codons: coordinates.append((i+1,0,codon))
00141             frame_coordinates.append(coordinates)
00142         return frame_coordinates

Here is the caller graph for this function:

def nextorf.NextOrf.Output (   self,
  CDS 
)

Definition at line 179 of file nextorf.py.

00179 
00180     def Output(self, CDS):
00181         out = self.options['output']
00182         seqs = (self.seq, self.rseq)
00183         n = len(self.seq)
00184         for start, stop, length, subs, strand in CDS:
00185             self.counter += 1
00186             if strand > 0: head = 'orf_%s:%s:%d:%d:%d' % (self.counter, self.header, strand, start,stop)
00187             if strand < 0: head = 'orf_%s:%s:%d:%d:%d' % (self.counter, self.header, strand, n-stop+1,n-start+1)
00188             if self.options['gc']:
00189                 head = '%s:%s' % (head, self.Gc2(subs.data))
00190                 
00191             if out == 'aa':
00192                 orf = subs.translate(table=self.genetic_code)
00193                 print self.ToFasta(head, orf.data)
00194             elif out == 'nt':
00195                 print self.ToFasta(head, subs.data)
00196             elif out == 'pos':
00197                 print head
00198                 
        

Here is the call graph for this function:

def nextorf.NextOrf.ReadFile (   self)

Definition at line 56 of file nextorf.py.

00056 
00057     def ReadFile(self):
00058         handle = open(self.file)
00059         for record in SeqIO.parse(handle, "fasta"):
00060             self.header = record.id
00061             frame_coordinates = ''
00062             dir = self.options['strand']
00063             plus = dir in ['both', 'plus']
00064             minus = dir in ['both', 'minus']
00065             start, stop = int(self.options['start']), int(self.options['stop'])
00066             s = str(record.seq).upper()
00067             if stop > 0:
00068                s = s[start:stop]
00069             else:
00070                s = s[start:]
00071             self.seq = Seq(s,IUPAC.ambiguous_dna)
00072             self.length = len(self.seq)
00073             self.rseq = None
00074             CDS = []
00075             if plus:
00076                 CDS.extend(self.GetCDS(self.seq))
00077             if minus:
00078                 self.rseq = self.seq.reverse_complement()
00079                 CDS.extend(self.GetCDS(self.rseq, strand = -1))
00080             self.Output(CDS)

Here is the call graph for this function:

def nextorf.NextOrf.ToFasta (   self,
  header,
  seq 
)

Definition at line 81 of file nextorf.py.

00081 
00082     def ToFasta(self, header, seq):
00083        seq = re.sub('(............................................................)','\\1\n',seq)
00084        return '>%s\n%s' % (header, seq)

Here is the caller graph for this function:


Member Data Documentation

Definition at line 53 of file nextorf.py.

Definition at line 50 of file nextorf.py.

Definition at line 51 of file nextorf.py.

Definition at line 59 of file nextorf.py.

Definition at line 71 of file nextorf.py.

Definition at line 49 of file nextorf.py.

Definition at line 72 of file nextorf.py.

Definition at line 70 of file nextorf.py.

Definition at line 52 of file nextorf.py.


The documentation for this class was generated from the following file: