<s> """<STR_LIT>""" <EOL> import re <EOL> import sys <EOL> import os <EOL> import os . path as op <EOL> import shutil <EOL> import logging <EOL> import string <EOL> from itertools import groupby , izip_longest <EOL> from Bio import SeqIO <EOL> from Bio . Seq import Seq <EOL> from Bio . SeqRecord import SeqRecord <EOL> from Bio . SeqUtils . CheckSum import seguid <EOL> from hashlib import md5 <EOL> from jcvi . formats . base import BaseFile , DictFile , must_open <EOL> from jcvi . formats . bed import Bed <EOL> from jcvi . utils . cbook import percentage <EOL> from jcvi . utils . table import write_csv <EOL> from jcvi . apps . console import red , green <EOL> from jcvi . apps . base import OptionParser , ActionDispatcher , need_update <EOL> class Fasta ( BaseFile , dict ) : <EOL> def __init__ ( self , filename , index = False , key_function = None , lazy = False ) : <EOL> super ( Fasta , self ) . __init__ ( filename ) <EOL> self . key_function = key_function <EOL> if lazy : <EOL> return <EOL> if index : <EOL> self . index = SeqIO . index ( filename , "<STR_LIT>" , <EOL> key_function = key_function ) <EOL> else : <EOL> _key_function = ( lambda rec : key_function ( rec . description ) ) if key_function else None <EOL> self . index = SeqIO . to_dict ( SeqIO . parse ( must_open ( filename ) , "<STR_LIT>" ) , <EOL> key_function = _key_function ) <EOL> def _key_function ( self , key ) : <EOL> return self . key_function ( key ) if self . key_function else key <EOL> def __len__ ( self ) : <EOL> return len ( self . index ) <EOL> def __contains__ ( self , key ) : <EOL> key = self . _key_function ( key ) <EOL> return key in self . index <EOL> def __getitem__ ( self , key ) : <EOL> key = self . _key_function ( key ) <EOL> rec = self . index [ key ] <EOL> return rec <EOL> def keys ( self ) : <EOL> return self . index . keys ( ) <EOL> def iterkeys ( self ) : <EOL> for k in self . index . iterkeys ( ) : <EOL> yield k <EOL> def iteritems ( self ) : <EOL> for k in self . iterkeys ( ) : <EOL> yield k , self [ k ] <EOL> def itersizes ( self ) : <EOL> for k in self . iterkeys ( ) : <EOL> yield k , len ( self [ k ] ) <EOL> def iteritems_ordered ( self ) : <EOL> for rec in SeqIO . parse ( must_open ( self . filename ) , "<STR_LIT>" ) : <EOL> yield rec . name , rec <EOL> def iterdescriptions_ordered ( self ) : <EOL> for k , rec in self . iteritems_ordered ( ) : <EOL> yield rec . description , rec <EOL> def iterkeys_ordered ( self ) : <EOL> for k , rec in self . iteritems_ordered ( ) : <EOL> yield k <EOL> def itersizes_ordered ( self ) : <EOL> for k , rec in self . iteritems_ordered ( ) : <EOL> yield k , len ( rec ) <EOL> def tostring ( self ) : <EOL> d = { } <EOL> for k , rec in self . iteritems ( ) : <EOL> d [ k ] = str ( rec . seq ) <EOL> return d <EOL> @ property <EOL> def totalsize ( self ) : <EOL> return sum ( size for k , size in self . itersizes ( ) ) <EOL> @ classmethod <EOL> def subseq ( cls , fasta , start = None , stop = None , strand = None ) : <EOL> """<STR_LIT>""" <EOL> start = start - <NUM_LIT:1> if start is not None else <NUM_LIT:0> <EOL> stop = stop if stop is not None else len ( fasta ) <EOL> if start < <NUM_LIT:0> : <EOL> msg = "<STR_LIT>" . format ( start + <NUM_LIT:1> , fasta . id ) <EOL> logging . error ( msg ) <EOL> start = <NUM_LIT:0> <EOL> if stop > len ( fasta ) : <EOL> msg = "<STR_LIT>" . format ( stop , fasta . id , len ( fasta ) ) <EOL> logging . error ( msg ) <EOL> stop = len ( fasta ) <EOL> seq = fasta . seq [ start : stop ] <EOL> if strand in ( - <NUM_LIT:1> , '<STR_LIT>' , '<STR_LIT:->' ) : <EOL> seq = seq . reverse_complement ( ) <EOL> return seq <EOL> def sequence ( self , f , asstring = True ) : <EOL> """<STR_LIT>""" <EOL> assert '<STR_LIT>' in f , "<STR_LIT>" <EOL> name = f [ '<STR_LIT>' ] <EOL> assert name in self , "<STR_LIT>" % ( f , self . filename ) <EOL> fasta = self [ f [ '<STR_LIT>' ] ] <EOL> seq = Fasta . subseq ( fasta , <EOL> f . get ( '<STR_LIT:start>' ) , f . get ( '<STR_LIT>' ) , f . get ( '<STR_LIT>' ) ) <EOL> if asstring : <EOL> return str ( seq ) <EOL> return seq <EOL> """<STR_LIT>""" <EOL> class ORFFinder ( object ) : <EOL> """<STR_LIT>""" <EOL> def __init__ ( self , seq , start = [ ] , stop = [ "<STR_LIT>" , "<STR_LIT>" , "<STR_LIT>" ] ) : <EOL> self . seq = seq . tostring ( ) . upper ( ) <EOL> self . start = start <EOL> self . stop = stop <EOL> self . result = [ "<STR_LIT:+>" , <NUM_LIT:0> , <NUM_LIT:0> , <NUM_LIT:0> , <NUM_LIT:0> ] <EOL> self . longest = <NUM_LIT:0> <EOL> self . size = len ( seq ) <EOL> def __str__ ( self ) : <EOL> strand , frame , start , end , length = self . result <EOL> start += <NUM_LIT:1> <EOL> if strand == '<STR_LIT:->' : <EOL> start , end = end , start <EOL> return "<STR_LIT>" . format ( start , end ) <EOL> @ property <EOL> def info ( self ) : <EOL> strand , frame , start , end , length = self . result <EOL> return "<STR_LIT:\t>" . join ( str ( x ) for x in ( strand , frame , start , end ) ) <EOL> def codons ( self , frame ) : <EOL> """<STR_LIT>""" <EOL> start = frame <EOL> while start + <NUM_LIT:3> <= self . size : <EOL> yield self . sequence [ start : start + <NUM_LIT:3> ] , start <EOL> start += <NUM_LIT:3> <EOL> def scan_sequence ( self , frame , direction ) : <EOL> """<STR_LIT>""" <EOL> orf_start = None <EOL> for c , index in self . codons ( frame ) : <EOL> if ( c not in self . stop and ( c in self . start or not self . start ) <EOL> and orf_start is None ) : <EOL> orf_start = index <EOL> elif c in self . stop and orf_start is not None : <EOL> self . _update_longest ( orf_start , index + <NUM_LIT:3> , direction , frame ) <EOL> orf_start = None <EOL> if orf_start is not None : <EOL> self . _update_longest ( orf_start , index + <NUM_LIT:3> , direction , frame ) <EOL> def _update_longest ( self , orf_start , index , direction , frame ) : <EOL> orf_end = index <EOL> L = orf_end - orf_start <EOL> if L > self . longest : <EOL> self . longest = L <EOL> self . result = [ direction , frame , orf_start , orf_end , L ] <EOL> def get_longest_orf ( self ) : <EOL> dirs = ( "<STR_LIT:+>" , "<STR_LIT:->" ) <EOL> for direction in dirs : <EOL> self . sequence = self . seq <EOL> if direction == "<STR_LIT:->" : <EOL> self . sequence = rc ( self . sequence ) <EOL> for frame in xrange ( <NUM_LIT:3> ) : <EOL> self . scan_sequence ( frame , direction ) <EOL> strand , frame , start , end , length = self . result <EOL> size = self . size <EOL> if strand == '<STR_LIT:->' : <EOL> start , end = size - end , size - start <EOL> self . result [ <NUM_LIT:2> : <NUM_LIT:4> ] = start , end <EOL> assert start <= end , self . result <EOL> if start == end : <EOL> return "<STR_LIT:N>" <EOL> orf = self . seq [ start : end ] <EOL> if strand == '<STR_LIT:->' : <EOL> orf = rc ( orf ) <EOL> assert len ( orf ) % <NUM_LIT:3> == <NUM_LIT:0> <EOL> return orf <EOL> class SequenceInfo ( object ) : <EOL> """<STR_LIT>""" <EOL> def __init__ ( self , filename , gapstats = False ) : <EOL> from jcvi . utils . cbook import SummaryStats <EOL> from jcvi . assembly . base import calculate_A50 <EOL> f = Fasta ( filename ) <EOL> self . filename = filename <EOL> self . header = "<STR_LIT>" . split ( "<STR_LIT:|>" ) <EOL> if gapstats : <EOL> self . header += [ "<STR_LIT>" ] <EOL> self . nseqs = len ( f ) <EOL> sizes = [ ] <EOL> gaps = [ ] <EOL> na = nc = ng = nt = <NUM_LIT:0> <EOL> for k , s in f . iteritems ( ) : <EOL> s = str ( s . seq ) . upper ( ) <EOL> sizes . append ( len ( s ) ) <EOL> na += s . count ( '<STR_LIT:A>' ) <EOL> nc += s . count ( '<STR_LIT:C>' ) <EOL> ng += s . count ( '<STR_LIT>' ) <EOL> nt += s . count ( '<STR_LIT:T>' ) <EOL> if gapstats : <EOL> gaps += list ( self . iter_gap_len ( s ) ) <EOL> self . real = real = na + nc + ng + nt <EOL> s = SummaryStats ( sizes ) <EOL> self . sum = s . sum <EOL> if gapstats : <EOL> self . gaps = len ( gaps ) <EOL> self . nn = self . sum - real <EOL> a50 , l50 , nn50 = calculate_A50 ( sizes ) <EOL> self . min = s . min <EOL> self . max = s . max <EOL> self . mean = int ( s . mean ) <EOL> self . n50 = l50 <EOL> self . data = [ self . filename , self . nseqs , <EOL> self . real , self . nn , self . sum , <EOL> self . min , self . max , self . n50 ] <EOL> if gapstats : <EOL> self . data += [ self . gaps ] <EOL> assert len ( self . header ) == len ( self . data ) <EOL> def iter_gap_len ( self , seq , mingap = <NUM_LIT:10> ) : <EOL> for gap , seq in groupby ( seq , lambda x : x == '<STR_LIT:N>' ) : <EOL> if not gap : <EOL> continue <EOL> gap_len = len ( list ( seq ) ) <EOL> if gap_len >= mingap : <EOL> yield len ( list ( seq ) ) <EOL> def rc ( s ) : <EOL> _complement = string . maketrans ( '<STR_LIT>' , '<STR_LIT>' ) <EOL> cs = s . translate ( _complement ) <EOL> return cs [ : : - <NUM_LIT:1> ] <EOL> def main ( ) : <EOL> actions = ( <EOL> ( '<STR_LIT>' , '<STR_LIT>' + '<STR_LIT>' ) , <EOL> ( '<STR_LIT>' , '<STR_LIT>' ) , <EOL> ( '<STR_LIT>' , '<STR_LIT>' ) , <EOL> ( '<STR_LIT:info>' , '<STR_LIT>' ) , <EOL> ( '<STR_LIT>' , "<STR_LIT>" ) , <EOL> ( '<STR_LIT>' , '<STR_LIT>' ) , <EOL> ( '<STR_LIT>' , '<STR_LIT>' ) , <EOL> ( '<STR_LIT>' , '<STR_LIT>' + '<STR_LIT>' ) , <EOL> ( '<STR_LIT>' , '<STR_LIT>' ) , <EOL> ( '<STR_LIT>' , '<STR_LIT>' ) , <EOL> ( '<STR_LIT>' , '<STR_LIT>' ) , <EOL> ( '<STR_LIT>' , '<STR_LIT>' ) , <EOL> ( '<STR_LIT>' , '<STR_LIT>' ) , <EOL> ( '<STR_LIT>' , '<STR_LIT>' ) , <EOL> ( '<STR_LIT>' , '<STR_LIT>' ) , <EOL> ( '<STR_LIT>' , '<STR_LIT>' ) , <EOL> ( '<STR_LIT>' , '<STR_LIT>' ) , <EOL> ( '<STR_LIT>' , '<STR_LIT>' + "<STR_LIT>" ) , <EOL> ( '<STR_LIT>' , '<STR_LIT>' ) , <EOL> ( '<STR_LIT>' , '<STR_LIT>' ) , <EOL> ( '<STR_LIT>' , '<STR_LIT>' ) , <EOL> ( '<STR_LIT>' , '<STR_LIT>' ) , <EOL> ( '<STR_LIT>' , '<STR_LIT>' ) , <EOL> ( '<STR_LIT>' , '<STR_LIT>' + '<STR_LIT>' ) , <EOL> ( '<STR_LIT>' , '<STR_LIT>' ) , <EOL> ( '<STR_LIT>' , '<STR_LIT>' ) , <EOL> ( '<STR_LIT>' , '<STR_LIT>' ) , <EOL> ( '<STR_LIT>' , '<STR_LIT>' ) , <EOL> ) <EOL> p = ActionDispatcher ( actions ) <EOL> p . dispatch ( globals ( ) ) <EOL> def trimsplit ( args ) : <EOL> """<STR_LIT>""" <EOL> from jcvi . utils . cbook import SummaryStats <EOL> p = OptionParser ( trimsplit . __doc__ ) <EOL> p . add_option ( "<STR_LIT>" , default = <NUM_LIT:1000> , type = "<STR_LIT:int>" , <EOL> help = "<STR_LIT>" ) <EOL> opts , args = p . parse_args ( args ) <EOL> if len ( args ) != <NUM_LIT:1> : <EOL> sys . exit ( not p . print_help ( ) ) <EOL> fastafile , = args <EOL> minlength = opts . minlength <EOL> fw = must_open ( fastafile . rsplit ( "<STR_LIT:.>" , <NUM_LIT:1> ) [ <NUM_LIT:0> ] + "<STR_LIT>" , "<STR_LIT:w>" ) <EOL> ntotal = <NUM_LIT:0> <EOL> removed = [ ] <EOL> Ns = [ ] <EOL> for name , seq in parse_fasta ( fastafile ) : <EOL> stretches = [ ] <EOL> ntotal += len ( seq ) <EOL> for lower , stretch in groupby ( seq , key = lambda x : x . islower ( ) ) : <EOL> stretch = "<STR_LIT>" . join ( stretch ) <EOL> if lower or len ( stretch ) < minlength : <EOL> removed . append ( len ( stretch ) ) <EOL> continue <EOL> for isN , s in groupby ( stretch , key = lambda x : x in "<STR_LIT>" ) : <EOL> s = "<STR_LIT>" . join ( s ) <EOL> if isN or len ( s ) < minlength : <EOL> Ns . append ( len ( s ) ) <EOL> continue <EOL> stretches . append ( s ) <EOL> for i , seq in enumerate ( stretches ) : <EOL> id = "<STR_LIT>" . format ( name . split ( "<STR_LIT:|>" ) [ <NUM_LIT:0> ] , i ) <EOL> s = SeqRecord ( Seq ( seq ) , id = id , description = "<STR_LIT>" ) <EOL> SeqIO . write ( [ s ] , fw , "<STR_LIT>" ) <EOL> fw . close ( ) <EOL> if removed : <EOL> logging . debug ( "<STR_LIT>" . format ( percentage ( sum ( removed ) , ntotal ) ) ) <EOL> print >> sys . stderr , SummaryStats ( removed ) <EOL> if Ns : <EOL> logging . debug ( "<STR_LIT>" . format ( percentage ( sum ( Ns ) , ntotal ) ) ) <EOL> print >> sys . stderr , SummaryStats ( Ns ) <EOL> def qual ( args ) : <EOL> """<STR_LIT>""" <EOL> from jcvi . formats . sizes import Sizes <EOL> p = OptionParser ( qual . __doc__ ) <EOL> p . add_option ( "<STR_LIT>" , default = <NUM_LIT> , type = "<STR_LIT:int>" , <EOL> help = "<STR_LIT>" ) <EOL> p . set_outfile ( ) <EOL> opts , args = p . parse_args ( args ) <EOL> if len ( args ) != <NUM_LIT:1> : <EOL> sys . exit ( not p . print_help ( ) ) <EOL> fastafile , = args <EOL> sizes = Sizes ( fastafile ) <EOL> qvchar = str ( opts . qv ) <EOL> fw = must_open ( opts . outfile , "<STR_LIT:w>" ) <EOL> total = <NUM_LIT:0> <EOL> for s , slen in sizes . iter_sizes ( ) : <EOL> print >> fw , "<STR_LIT:>>" + s <EOL> print >> fw , "<STR_LIT:U+0020>" . join ( [ qvchar ] * slen ) <EOL> total += <NUM_LIT:1> <EOL> fw . close ( ) <EOL> logging . debug ( "<STR_LIT>" . format ( total , opts . outfile ) ) <EOL> def info ( args ) : <EOL> """<STR_LIT>""" <EOL> p = OptionParser ( info . __doc__ ) <EOL> p . add_option ( "<STR_LIT>" , default = False , action = "<STR_LIT:store_true>" , <EOL> help = "<STR_LIT>" ) <EOL> p . set_table ( ) <EOL> p . set_outfile ( ) <EOL> opts , args = p . parse_args ( args ) <EOL> if len ( args ) == <NUM_LIT:0> : <EOL> sys . exit ( not p . print_help ( ) ) <EOL> fastafiles = args <EOL> data = [ ] <EOL> for f in fastafiles : <EOL> s = SequenceInfo ( f , gapstats = opts . gaps ) <EOL> data . append ( s . data ) <EOL> write_csv ( s . header , data , sep = opts . sep , <EOL> filename = opts . outfile , align = opts . align ) <EOL> def fromtab ( args ) : <EOL> """<STR_LIT>""" <EOL> p = OptionParser ( fromtab . __doc__ ) <EOL> p . set_sep ( sep = None ) <EOL> p . add_option ( "<STR_LIT>" , default = False , action = "<STR_LIT:store_true>" , <EOL> help = "<STR_LIT>" ) <EOL> p . add_option ( "<STR_LIT>" , <EOL> help = "<STR_LIT>" ) <EOL> opts , args = p . parse_args ( args ) <EOL> if len ( args ) != <NUM_LIT:2> : <EOL> sys . exit ( not p . print_help ( ) ) <EOL> tabfile , fastafile = args <EOL> sep = opts . sep <EOL> replace = opts . replace <EOL> fp = must_open ( tabfile ) <EOL> fw = must_open ( fastafile , "<STR_LIT:w>" ) <EOL> nseq = <NUM_LIT:0> <EOL> if opts . noheader : <EOL> fp . next ( ) <EOL> for row in fp : <EOL> row = row . strip ( ) <EOL> if not row or row [ <NUM_LIT:0> ] == '<STR_LIT:#>' : <EOL> continue <EOL> name , seq = row . rsplit ( sep , <NUM_LIT:1> ) <EOL> if replace : <EOL> name = name . replace ( "<STR_LIT:U+0020>" , replace ) <EOL> print >> fw , "<STR_LIT>" . format ( name , seq ) <EOL> nseq += <NUM_LIT:1> <EOL> fw . close ( ) <EOL> logging . debug ( "<STR_LIT>" . format ( nseq , fastafile ) ) <EOL> def longestorf ( args ) : <EOL> """<STR_LIT>""" <EOL> p = OptionParser ( longestorf . __doc__ ) <EOL> p . add_option ( "<STR_LIT>" , action = "<STR_LIT:store_true>" , <EOL> help = "<STR_LIT>" ) <EOL> opts , args = p . parse_args ( args ) <EOL> if len ( args ) != <NUM_LIT:1> : <EOL> sys . exit ( not p . print_help ( ) ) <EOL> fastafile , = args <EOL> pf = fastafile . rsplit ( "<STR_LIT:.>" , <NUM_LIT:1> ) [ <NUM_LIT:0> ] <EOL> orffile = pf + "<STR_LIT>" <EOL> idsfile = None <EOL> if opts . ids : <EOL> idsfile = pf + "<STR_LIT>" <EOL> fwids = open ( idsfile , "<STR_LIT:w>" ) <EOL> f = Fasta ( fastafile , lazy = True ) <EOL> fw = must_open ( orffile , "<STR_LIT:w>" ) <EOL> before , after = <NUM_LIT:0> , <NUM_LIT:0> <EOL> for name , rec in f . iteritems_ordered ( ) : <EOL> cds = rec . seq <EOL> before += len ( cds ) <EOL> orf = ORFFinder ( cds ) <EOL> lorf = orf . get_longest_orf ( ) <EOL> newcds = Seq ( lorf ) <EOL> after += len ( newcds ) <EOL> newrec = SeqRecord ( newcds , id = name , description = rec . description ) <EOL> SeqIO . write ( [ newrec ] , fw , "<STR_LIT>" ) <EOL> if idsfile : <EOL> print >> fwids , "<STR_LIT:\t>" . join ( ( name , orf . info ) ) <EOL> fw . close ( ) <EOL> if idsfile : <EOL> fwids . close ( ) <EOL> logging . debug ( "<STR_LIT>" . format ( orffile , percentage ( after , before ) ) ) <EOL> return orffile <EOL> def ispcr ( args ) : <EOL> """<STR_LIT>""" <EOL> from jcvi . utils . iter import grouper <EOL> p = OptionParser ( ispcr . __doc__ ) <EOL> p . add_option ( "<STR_LIT>" , dest = "<STR_LIT>" , default = <NUM_LIT:1> , type = "<STR_LIT:int>" , <EOL> help = "<STR_LIT>" ) <EOL> opts , args = p . parse_args ( args ) <EOL> if len ( args ) != <NUM_LIT:1> : <EOL> sys . exit ( not p . print_help ( ) ) <EOL> fastafile , = args <EOL> ispcrfile = fastafile + "<STR_LIT>" <EOL> fw = open ( ispcrfile , "<STR_LIT:w>" ) <EOL> N = opts . rclip <EOL> strip_name = lambda x : x [ : - N ] if N else str <EOL> npairs = <NUM_LIT:0> <EOL> fastaiter = SeqIO . parse ( fastafile , "<STR_LIT>" ) <EOL> for a , b in grouper ( fastaiter , <NUM_LIT:2> ) : <EOL> aid , bid = [ strip_name ( x ) for x in ( a . id , b . id ) ] <EOL> assert aid == bid , "<STR_LIT>" . format ( ( aid , bid ) ) <EOL> print >> fw , "<STR_LIT:\t>" . join ( ( aid , str ( a . seq ) , str ( b . seq ) ) ) <EOL> npairs += <NUM_LIT:1> <EOL> fw . close ( ) <EOL> logging . debug ( "<STR_LIT>" . format ( npairs , ispcrfile ) ) <EOL> def parse_fasta ( infile , upper = False ) : <EOL> '''<STR_LIT>''' <EOL> try : <EOL> fp = must_open ( infile ) <EOL> except : <EOL> fp = infile <EOL> fa_iter = ( x [ <NUM_LIT:1> ] for x in groupby ( fp , lambda row : row [ <NUM_LIT:0> ] == '<STR_LIT:>>' ) ) <EOL> for header in fa_iter : <EOL> header = header . next ( ) <EOL> if header [ <NUM_LIT:0> ] != '<STR_LIT:>>' : <EOL> continue <EOL> header = header . strip ( ) [ <NUM_LIT:1> : ] <EOL> seq = "<STR_LIT>" . join ( s . strip ( ) for s in fa_iter . next ( ) ) <EOL> if upper : <EOL> seq = seq . upper ( ) <EOL> yield header , seq <EOL> def iter_clean_fasta ( fastafile ) : <EOL> for header , seq in parse_fasta ( fastafile ) : <EOL> seq = "<STR_LIT>" . join ( x for x in seq if x in string . letters or x == '<STR_LIT:*>' ) <EOL> yield header , seq <EOL> def iter_canonical_fasta ( fastafile ) : <EOL> canonical = "<STR_LIT>" <EOL> totalbad = <NUM_LIT:0> <EOL> for header , seq in parse_fasta ( fastafile ) : <EOL> badcounts = sum ( <NUM_LIT:1> for x in seq if x not in canonical ) <EOL> seq = "<STR_LIT>" . join ( ( x if x in canonical else '<STR_LIT:N>' ) for x in seq ) <EOL> totalbad += badcounts <EOL> yield header , seq <EOL> logging . debug ( "<STR_LIT>" . format ( totalbad ) ) <EOL> def fancyprint ( fw , seq , width = <NUM_LIT> , chunk = <NUM_LIT:10> ) : <EOL> from jcvi . utils . iter import grouper <EOL> assert width % chunk == <NUM_LIT:0> <EOL> nchunks = width / chunk <EOL> seqlen = len ( seq ) <EOL> maxchar = len ( str ( seqlen ) ) <EOL> s = [ "<STR_LIT>" . join ( x ) for x in grouper ( seq , chunk , fillvalue = "<STR_LIT>" ) ] <EOL> s = [ "<STR_LIT:U+0020>" . join ( x ) for x in grouper ( s , nchunks , fillvalue = "<STR_LIT>" ) ] <EOL> for a , b in zip ( range ( <NUM_LIT:1> , len ( seq ) , width ) , s ) : <EOL> b = b . rstrip ( ) <EOL> a = str ( a ) . rjust ( maxchar , "<STR_LIT:U+0020>" ) <EOL> print >> fw , "<STR_LIT:U+0020>" . join ( ( a , b ) ) <EOL> def clean ( args ) : <EOL> """<STR_LIT>""" <EOL> p = OptionParser ( clean . __doc__ ) <EOL> p . add_option ( "<STR_LIT>" , default = False , action = "<STR_LIT:store_true>" , <EOL> help = "<STR_LIT>" ) <EOL> p . add_option ( "<STR_LIT>" , default = False , action = "<STR_LIT:store_true>" , <EOL> help = "<STR_LIT>" ) <EOL> p . set_outfile ( ) <EOL> opts , args = p . parse_args ( args ) <EOL> if len ( args ) != <NUM_LIT:1> : <EOL> sys . exit ( not p . print_help ( ) ) <EOL> fastafile , = args <EOL> fw = must_open ( opts . outfile , "<STR_LIT:w>" ) <EOL> if opts . fancy : <EOL> for header , seq in iter_clean_fasta ( fastafile ) : <EOL> print >> fw , "<STR_LIT:>>" + header <EOL> fancyprint ( fw , seq ) <EOL> return <NUM_LIT:0> <EOL> iterator = iter_canonical_fasta if opts . canonical else iter_clean_fasta <EOL> for header , seq in iterator ( fastafile ) : <EOL> seq = Seq ( seq ) <EOL> s = SeqRecord ( seq , id = header , description = "<STR_LIT>" ) <EOL> SeqIO . write ( [ s ] , fw , "<STR_LIT>" ) <EOL> def translate ( args ) : <EOL> """<STR_LIT>""" <EOL> transl_tables = [ str ( x ) for x in xrange ( <NUM_LIT:1> , <NUM_LIT> ) ] <EOL> p = OptionParser ( translate . __doc__ ) <EOL> p . add_option ( "<STR_LIT>" , default = False , action = "<STR_LIT:store_true>" , <EOL> help = "<STR_LIT>" <EOL> "<STR_LIT>" ) <EOL> p . add_option ( "<STR_LIT>" , default = False , action = "<STR_LIT:store_true>" , <EOL> help = "<STR_LIT>" ) <EOL> p . add_option ( "<STR_LIT>" , default = <NUM_LIT:1> , choices = transl_tables , <EOL> help = "<STR_LIT>" ) <EOL> p . set_outfile ( ) <EOL> opts , args = p . parse_args ( args ) <EOL> if len ( args ) != <NUM_LIT:1> : <EOL> sys . exit ( not p . print_help ( ) ) <EOL> cdsfasta , = args <EOL> if opts . longest : <EOL> cdsfasta = longestorf ( [ cdsfasta ] ) <EOL> f = Fasta ( cdsfasta , lazy = True ) <EOL> outfile = opts . outfile <EOL> fw = must_open ( outfile , "<STR_LIT:w>" ) <EOL> if opts . ids : <EOL> idsfile = cdsfasta . rsplit ( "<STR_LIT:.>" , <NUM_LIT:1> ) [ <NUM_LIT:0> ] + "<STR_LIT>" <EOL> ids = open ( idsfile , "<STR_LIT:w>" ) <EOL> else : <EOL> ids = None <EOL> five_prime_missing = three_prime_missing = <NUM_LIT:0> <EOL> contain_ns = complete = cannot_translate = total = <NUM_LIT:0> <EOL> for name , rec in f . iteritems_ordered ( ) : <EOL> cds = rec . seq <EOL> cdslen = len ( cds ) <EOL> peplen = cdslen / <NUM_LIT:3> <EOL> total += <NUM_LIT:1> <EOL> pep = "<STR_LIT>" <EOL> for i in xrange ( <NUM_LIT:3> ) : <EOL> newcds = cds [ i : i + peplen * <NUM_LIT:3> ] <EOL> newpep = newcds . translate ( table = opts . table ) <EOL> if len ( newpep . split ( "<STR_LIT:*>" ) [ <NUM_LIT:0> ] ) > len ( pep . split ( "<STR_LIT:*>" ) [ <NUM_LIT:0> ] ) : <EOL> pep = newpep <EOL> labels = [ ] <EOL> if "<STR_LIT:*>" in pep . rstrip ( "<STR_LIT:*>" ) : <EOL> logging . error ( "<STR_LIT>" . format ( name ) ) <EOL> cannot_translate += <NUM_LIT:1> <EOL> labels . append ( "<STR_LIT>" ) <EOL> contains_start = pep . startswith ( "<STR_LIT:M>" ) <EOL> contains_stop = pep . endswith ( "<STR_LIT:*>" ) <EOL> contains_ns = "<STR_LIT:X>" in pep <EOL> start_ns = pep . startswith ( "<STR_LIT:X>" ) <EOL> end_ns = pep . endswith ( "<STR_LIT:X>" ) <EOL> if not contains_start : <EOL> five_prime_missing += <NUM_LIT:1> <EOL> labels . append ( "<STR_LIT>" ) <EOL> if not contains_stop : <EOL> three_prime_missing += <NUM_LIT:1> <EOL> labels . append ( "<STR_LIT>" ) <EOL> if contains_ns : <EOL> contain_ns += <NUM_LIT:1> <EOL> labels . append ( "<STR_LIT>" ) <EOL> if contains_start and contains_stop : <EOL> complete += <NUM_LIT:1> <EOL> labels . append ( "<STR_LIT>" ) <EOL> if start_ns : <EOL> labels . append ( "<STR_LIT>" ) <EOL> if end_ns : <EOL> labels . append ( "<STR_LIT>" ) <EOL> if ids : <EOL> print >> ids , "<STR_LIT:\t>" . join ( ( name , "<STR_LIT:U+002C>" . join ( labels ) ) ) <EOL> peprec = SeqRecord ( pep , id = name , description = rec . description ) <EOL> SeqIO . write ( [ peprec ] , fw , "<STR_LIT>" ) <EOL> fw . flush ( ) <EOL> print >> sys . stderr , "<STR_LIT>" . format ( percentage ( complete , total ) ) <EOL> print >> sys . stderr , "<STR_LIT>" . format ( percentage ( five_prime_missing , total ) ) <EOL> print >> sys . stderr , "<STR_LIT>" . format ( percentage ( three_prime_missing , total ) ) <EOL> print >> sys . stderr , "<STR_LIT>" . format ( percentage ( contain_ns , total ) ) <EOL> if cannot_translate : <EOL> print >> sys . stderr , "<STR_LIT>" . format ( percentage ( cannot_translate , total ) ) <EOL> fw . close ( ) <EOL> return cdsfasta , outfile <EOL> def filter ( args ) : <EOL> """<STR_LIT>""" <EOL> p = OptionParser ( filter . __doc__ ) <EOL> p . add_option ( "<STR_LIT>" , default = False , action = "<STR_LIT:store_true>" , <EOL> help = "<STR_LIT>" ) <EOL> p . set_outfile ( ) <EOL> opts , args = p . parse_args ( args ) <EOL> if len ( args ) != <NUM_LIT:2> : <EOL> sys . exit ( not p . print_help ( ) ) <EOL> fastafile , cutoff = args <EOL> try : <EOL> cutoff = int ( cutoff ) <EOL> except ValueError : <EOL> sys . exit ( not p . print_help ( ) ) <EOL> f = Fasta ( fastafile , lazy = True ) <EOL> fw = must_open ( opts . outfile , "<STR_LIT:w>" ) <EOL> for name , rec in f . iteritems_ordered ( ) : <EOL> if opts . less and len ( rec ) >= cutoff : <EOL> continue <EOL> if ( not opts . less ) and len ( rec ) < cutoff : <EOL> continue <EOL> SeqIO . write ( [ rec ] , fw , "<STR_LIT>" ) <EOL> fw . flush ( ) <EOL> return fw . name <EOL> def pool ( args ) : <EOL> """<STR_LIT>""" <EOL> from jcvi . formats . base import longest_unique_prefix <EOL> p = OptionParser ( pool . __doc__ ) <EOL> opts , args = p . parse_args ( args ) <EOL> if len ( args ) < <NUM_LIT:1> : <EOL> sys . exit ( not p . print_help ( ) ) <EOL> for fastafile in args : <EOL> pf = longest_unique_prefix ( fastafile , args ) <EOL> print >> sys . stderr , fastafile , "<STR_LIT>" , pf <EOL> prefixopt = "<STR_LIT>" . format ( pf ) <EOL> format ( [ fastafile , "<STR_LIT>" , prefixopt ] ) <EOL> def ids ( args ) : <EOL> """<STR_LIT>""" <EOL> p = OptionParser ( ids . __doc__ ) <EOL> p . add_option ( "<STR_LIT>" , default = None , <EOL> help = "<STR_LIT>" ) <EOL> p . add_option ( "<STR_LIT>" , default = False , action = "<STR_LIT:store_true>" , <EOL> help = "<STR_LIT>" ) <EOL> p . set_outfile ( ) <EOL> opts , args = p . parse_args ( args ) <EOL> if len ( args ) < <NUM_LIT:1> : <EOL> sys . exit ( not p . print_help ( ) ) <EOL> until = opts . until <EOL> fw = must_open ( opts . outfile , "<STR_LIT:w>" ) <EOL> for row in must_open ( args ) : <EOL> if row [ <NUM_LIT:0> ] == "<STR_LIT:>>" : <EOL> row = row [ <NUM_LIT:1> : ] . rstrip ( ) <EOL> if until : <EOL> row = row . split ( until ) [ <NUM_LIT:0> ] <EOL> atoms = row . split ( None , <NUM_LIT:1> ) <EOL> if opts . description : <EOL> outrow = "<STR_LIT:\t>" . join ( atoms ) <EOL> else : <EOL> outrow = atoms [ <NUM_LIT:0> ] <EOL> print >> fw , outrow <EOL> fw . close ( ) <EOL> def sort ( args ) : <EOL> """<STR_LIT>""" <EOL> p = OptionParser ( sort . __doc__ ) <EOL> p . add_option ( "<STR_LIT>" , default = False , action = "<STR_LIT:store_true>" , <EOL> help = "<STR_LIT>" ) <EOL> opts , args = p . parse_args ( args ) <EOL> if len ( args ) != <NUM_LIT:1> : <EOL> sys . exit ( p . print_help ( ) ) <EOL> fastafile , = args <EOL> sortedfastafile = fastafile . rsplit ( "<STR_LIT:.>" , <NUM_LIT:1> ) [ <NUM_LIT:0> ] + "<STR_LIT>" <EOL> f = Fasta ( fastafile , index = False ) <EOL> fw = must_open ( sortedfastafile , "<STR_LIT:w>" ) <EOL> if opts . sizes : <EOL> sortlist = sorted ( f . itersizes ( ) , key = lambda x : ( - x [ <NUM_LIT:1> ] , x [ <NUM_LIT:0> ] ) ) <EOL> logging . debug ( "<STR_LIT>" . format ( sortlist [ <NUM_LIT:0> ] , sortlist [ - <NUM_LIT:1> ] ) ) <EOL> sortlist = [ x for x , s in sortlist ] <EOL> else : <EOL> sortlist = sorted ( f . iterkeys ( ) ) <EOL> for key in sortlist : <EOL> rec = f [ key ] <EOL> SeqIO . write ( [ rec ] , fw , "<STR_LIT>" ) <EOL> logging . debug ( "<STR_LIT>" . format ( sortedfastafile ) ) <EOL> fw . close ( ) <EOL> return sortedfastafile <EOL> def join ( args ) : <EOL> """<STR_LIT>""" <EOL> from jcvi . formats . agp import OO , Phases , build <EOL> from jcvi . formats . sizes import Sizes <EOL> p = OptionParser ( join . __doc__ ) <EOL> p . add_option ( "<STR_LIT>" , default = None , <EOL> help = "<STR_LIT>" ) <EOL> p . add_option ( "<STR_LIT>" , default = <NUM_LIT:100> , type = "<STR_LIT:int>" , <EOL> help = "<STR_LIT>" ) <EOL> p . add_option ( "<STR_LIT>" , default = "<STR_LIT>" , <EOL> help = "<STR_LIT>" ) <EOL> p . add_option ( "<STR_LIT>" , default = "<STR_LIT>" , <EOL> help = "<STR_LIT>" ) <EOL> p . add_option ( "<STR_LIT>" , help = "<STR_LIT>" ) <EOL> opts , args = p . parse_args ( args ) <EOL> nargs = len ( args ) <EOL> if nargs not in ( <NUM_LIT:1> , <NUM_LIT:2> ) : <EOL> sys . exit ( not p . print_help ( ) ) <EOL> if nargs == <NUM_LIT:2> : <EOL> fastafile , phasefile = args <EOL> phases = DictFile ( phasefile ) <EOL> phases = dict ( ( a , Phases [ int ( b ) ] ) for a , b in phases . items ( ) ) <EOL> else : <EOL> fastafile , = args <EOL> phases = { } <EOL> sizes = Sizes ( fastafile ) <EOL> prefix = fastafile . rsplit ( "<STR_LIT:.>" , <NUM_LIT:1> ) [ <NUM_LIT:0> ] <EOL> agpfile = prefix + "<STR_LIT>" <EOL> newid = opts . newid <EOL> oo = opts . oo <EOL> o = OO ( oo , sizes . mapping ) <EOL> if oo : <EOL> seen = o . contigs <EOL> logging . debug ( "<STR_LIT>" . format ( len ( sizes ) , len ( seen ) , oo ) ) <EOL> for ctg , size in sizes . iter_sizes ( ) : <EOL> if ctg in seen : <EOL> continue <EOL> o . add ( ctg , ctg , size ) <EOL> else : <EOL> if newid : <EOL> for ctg , size in sizes . iter_sizes ( ) : <EOL> o . add ( newid , ctg , size ) <EOL> else : <EOL> for scaffold_number , ( ctg , size ) in enumerate ( sizes . iter_sizes ( ) ) : <EOL> object_id = "<STR_LIT>" . format ( scaffold_number + <NUM_LIT:1> ) <EOL> o . add ( object_id , ctg , size ) <EOL> fw = open ( agpfile , "<STR_LIT:w>" ) <EOL> o . write_AGP ( fw , gapsize = opts . gapsize , gaptype = opts . gaptype , <EOL> evidence = opts . evidence , phases = phases ) <EOL> fw . close ( ) <EOL> joinedfastafile = prefix + "<STR_LIT>" <EOL> build ( [ agpfile , fastafile , joinedfastafile ] ) <EOL> return joinedfastafile <EOL> def summary ( args ) : <EOL> """<STR_LIT>""" <EOL> from jcvi . utils . natsort import natsort_key <EOL> p = OptionParser ( summary . __doc__ ) <EOL> p . add_option ( "<STR_LIT>" , default = "<STR_LIT>" , <EOL> help = "<STR_LIT>" ) <EOL> p . add_option ( "<STR_LIT>" , <EOL> help = "<STR_LIT>" ) <EOL> p . set_outfile ( ) <EOL> opts , args = p . parse_args ( args ) <EOL> if len ( args ) == <NUM_LIT:0> : <EOL> sys . exit ( not p . print_help ( ) ) <EOL> idsfile = opts . ids <EOL> header = "<STR_LIT>" . split ( ) <EOL> if idsfile : <EOL> idsfile = open ( idsfile , "<STR_LIT:w>" ) <EOL> nids = <NUM_LIT:0> <EOL> data = [ ] <EOL> for fastafile in args : <EOL> for rec in SeqIO . parse ( fastafile , "<STR_LIT>" ) : <EOL> seqlen = len ( rec ) <EOL> nns = rec . seq . count ( '<STR_LIT:n>' ) + rec . seq . count ( '<STR_LIT:N>' ) <EOL> reals = seqlen - nns <EOL> pct = reals * <NUM_LIT> / seqlen <EOL> pctreal = "<STR_LIT>" . format ( pct ) <EOL> if idsfile and pct < <NUM_LIT:50> : <EOL> nids += <NUM_LIT:1> <EOL> print >> idsfile , rec . id <EOL> data . append ( ( rec . id , reals , nns , seqlen , pctreal ) ) <EOL> data . sort ( key = natsort_key ) <EOL> ids , reals , nns , seqlen , pctreal = zip ( * data ) <EOL> reals = sum ( reals ) <EOL> nns = sum ( nns ) <EOL> seqlen = sum ( seqlen ) <EOL> pctreal = "<STR_LIT>" . format ( reals * <NUM_LIT> / seqlen ) <EOL> data . append ( ( "<STR_LIT>" , reals , nns , seqlen , pctreal ) ) <EOL> write_csv ( header , data , sep = "<STR_LIT:U+0020>" , filename = opts . outfile , thousands = True ) <EOL> if idsfile : <EOL> logging . debug ( "<STR_LIT>" . format ( nids , idsfile . name ) ) <EOL> idsfile . close ( ) <EOL> return reals , nns , seqlen <EOL> def format ( args ) : <EOL> """<STR_LIT>""" <EOL> sequential_choices = ( "<STR_LIT:replace>" , "<STR_LIT>" , "<STR_LIT>" ) <EOL> p = OptionParser ( format . __doc__ ) <EOL> p . add_option ( "<STR_LIT>" , default = False , action = "<STR_LIT:store_true>" , <EOL> help = "<STR_LIT>" ) <EOL> p . add_option ( "<STR_LIT>" , default = None , choices = sequential_choices , <EOL> help = "<STR_LIT>" ) <EOL> p . add_option ( "<STR_LIT>" , default = <NUM_LIT:0> , type = "<STR_LIT:int>" , <EOL> help = "<STR_LIT>" ) <EOL> p . add_option ( "<STR_LIT>" , default = <NUM_LIT:0> , type = "<STR_LIT:int>" , <EOL> help = "<STR_LIT>" ) <EOL> p . add_option ( "<STR_LIT>" , default = False , action = "<STR_LIT:store_true>" , <EOL> help = "<STR_LIT>" ) <EOL> p . add_option ( "<STR_LIT>" , default = None , <EOL> help = "<STR_LIT>" ) <EOL> p . add_option ( "<STR_LIT>" , default = <NUM_LIT:0> , type = "<STR_LIT:int>" , <EOL> help = "<STR_LIT>" ) <EOL> p . add_option ( "<STR_LIT>" , default = False , action = "<STR_LIT:store_true>" , <EOL> help = "<STR_LIT>" ) <EOL> p . add_option ( "<STR_LIT>" , help = "<STR_LIT>" ) <EOL> p . add_option ( "<STR_LIT>" , help = "<STR_LIT>" ) <EOL> p . add_option ( "<STR_LIT>" , default = False , action = "<STR_LIT:store_true>" , <EOL> help = "<STR_LIT>" ) <EOL> p . add_option ( "<STR_LIT>" , help = "<STR_LIT>" ) <EOL> p . add_option ( "<STR_LIT>" , help = "<STR_LIT>" <EOL> "<STR_LIT>" ) <EOL> p . add_option ( "<STR_LIT>" , help = "<STR_LIT>" ) <EOL> p . add_option ( "<STR_LIT>" , default = False , action = "<STR_LIT:store_true>" , <EOL> help = "<STR_LIT>" ) <EOL> p . add_option ( "<STR_LIT>" , default = False , action = "<STR_LIT:store_true>" , <EOL> help = "<STR_LIT>" ) <EOL> opts , args = p . parse_args ( args ) <EOL> if len ( args ) != <NUM_LIT:2> : <EOL> sys . exit ( not p . print_help ( ) ) <EOL> infasta , outfasta = args <EOL> gb = opts . gb <EOL> pairs = opts . pairs <EOL> prefix = opts . prefix <EOL> suffix = opts . suffix <EOL> noversion = opts . noversion <EOL> sequential = opts . sequential <EOL> sequentialoffset = opts . sequentialoffset <EOL> sep = opts . sep <EOL> idx = opts . index <EOL> mapfile = opts . switch <EOL> annotfile = opts . annotation <EOL> desc = not opts . nodesc <EOL> idsfile = opts . ids <EOL> idsfile = open ( idsfile , "<STR_LIT:w>" ) if idsfile else None <EOL> upper = opts . upper <EOL> if mapfile : <EOL> mapping = DictFile ( mapfile , delimiter = "<STR_LIT:\t>" ) <EOL> if annotfile : <EOL> annotation = DictFile ( annotfile , delimiter = "<STR_LIT:\t>" ) <EOL> fp = SeqIO . parse ( must_open ( infasta ) , "<STR_LIT>" ) <EOL> fw = must_open ( outfasta , "<STR_LIT:w>" ) <EOL> for i , rec in enumerate ( fp ) : <EOL> origid = rec . id <EOL> description = rec . description . replace ( origid , "<STR_LIT>" ) . strip ( ) <EOL> if sep : <EOL> rec . id = rec . description . split ( sep ) [ idx ] . strip ( ) <EOL> if gb : <EOL> atoms = rec . id . split ( "<STR_LIT:|>" ) <EOL> if len ( atoms ) >= <NUM_LIT:3> : <EOL> rec . id = atoms [ <NUM_LIT:3> ] <EOL> elif len ( atoms ) == <NUM_LIT:2> : <EOL> rec . id = atoms [ <NUM_LIT:1> ] <EOL> if pairs : <EOL> id = "<STR_LIT>" if ( i % <NUM_LIT:2> == <NUM_LIT:0> ) else "<STR_LIT>" <EOL> rec . id += id <EOL> if noversion : <EOL> rec . id = rec . id . rsplit ( "<STR_LIT:.>" , <NUM_LIT:1> ) [ <NUM_LIT:0> ] <EOL> if sequential : <EOL> rec . id = "<STR_LIT>" . format ( sequentialoffset , opts . pad0 ) <EOL> if sequential == "<STR_LIT>" : <EOL> rec . id = "<STR_LIT>" . format ( rec . id , origid ) <EOL> elif sequential == "<STR_LIT>" : <EOL> rec . id = "<STR_LIT>" . format ( origid , rec . id ) <EOL> sequentialoffset += <NUM_LIT:1> <EOL> if opts . template : <EOL> template , dir , lib = [ x . split ( "<STR_LIT:=>" ) [ - <NUM_LIT:1> ] for x in <EOL> rec . description . split ( ) [ <NUM_LIT:1> : <NUM_LIT:4> ] ] <EOL> rec . id = "<STR_LIT>" . format ( lib , template , dir ) <EOL> if mapfile : <EOL> if origid in mapping : <EOL> rec . id = mapping [ origid ] <EOL> else : <EOL> logging . error ( "<STR_LIT>" . format ( origid , mapfile ) ) <EOL> if prefix : <EOL> rec . id = prefix + rec . id <EOL> if suffix : <EOL> rec . id += suffix <EOL> if annotfile : <EOL> rec . description = annotation . get ( origid , "<STR_LIT>" ) if not mapfile else annotation . get ( rec . id , "<STR_LIT>" ) <EOL> else : <EOL> rec . description = description if desc else "<STR_LIT>" <EOL> if idsfile : <EOL> print >> idsfile , "<STR_LIT:\t>" . join ( ( origid , rec . id ) ) <EOL> if upper : <EOL> rec . seq = rec . seq . upper ( ) <EOL> SeqIO . write ( rec , fw , "<STR_LIT>" ) <EOL> if idsfile : <EOL> logging . debug ( "<STR_LIT>" . format ( idsfile . name ) ) <EOL> idsfile . close ( ) <EOL> def print_first_difference ( arec , brec , ignore_case = False , ignore_N = False , <EOL> rc = False , report_match = True ) : <EOL> """<STR_LIT>""" <EOL> plus_match = _print_first_difference ( arec , brec , ignore_case = ignore_case , <EOL> ignore_N = ignore_N , report_match = report_match ) <EOL> if rc and not plus_match : <EOL> logging . debug ( "<STR_LIT>" % brec . id ) <EOL> brec . seq = brec . seq . reverse_complement ( ) <EOL> minus_match = _print_first_difference ( arec , brec , <EOL> ignore_case = ignore_case , ignore_N = ignore_N , <EOL> report_match = report_match ) <EOL> return minus_match <EOL> else : <EOL> return plus_match <EOL> def _print_first_difference ( arec , brec , ignore_case = False , ignore_N = False , <EOL> report_match = True ) : <EOL> """<STR_LIT>""" <EOL> aseq , bseq = arec . seq , brec . seq <EOL> asize , bsize = len ( aseq ) , len ( bseq ) <EOL> matched = True <EOL> for i , ( a , b ) in enumerate ( izip_longest ( aseq , bseq ) ) : <EOL> if ignore_case and None not in ( a , b ) : <EOL> a , b = a . upper ( ) , b . upper ( ) <EOL> if ignore_N and ( '<STR_LIT:N>' in ( a , b ) or '<STR_LIT:X>' in ( a , b ) ) : <EOL> continue <EOL> if a != b : <EOL> matched = False <EOL> break <EOL> if i + <NUM_LIT:1> == asize and matched : <EOL> if report_match : <EOL> print green ( "<STR_LIT>" ) <EOL> match = True <EOL> else : <EOL> print red ( "<STR_LIT>" ) <EOL> snippet_size = <NUM_LIT:20> <EOL> print red ( "<STR_LIT>" % ( i + <NUM_LIT:1> ) ) <EOL> begin = max ( i - snippet_size , <NUM_LIT:0> ) <EOL> aend = min ( i + snippet_size , asize ) <EOL> bend = min ( i + snippet_size , bsize ) <EOL> print red ( aseq [ begin : i ] + "<STR_LIT:|>" + aseq [ i : aend ] ) <EOL> print red ( bseq [ begin : i ] + "<STR_LIT:|>" + bseq [ i : bend ] ) <EOL> match = False <EOL> return match <EOL> def diff ( args ) : <EOL> """<STR_LIT>""" <EOL> from jcvi . utils . table import banner <EOL> p = OptionParser ( diff . __doc__ ) <EOL> p . add_option ( "<STR_LIT>" , default = False , action = "<STR_LIT:store_true>" , <EOL> help = "<STR_LIT>" ) <EOL> p . add_option ( "<STR_LIT>" , default = False , action = "<STR_LIT:store_true>" , <EOL> help = "<STR_LIT>" ) <EOL> p . add_option ( "<STR_LIT>" , default = False , action = "<STR_LIT:store_true>" , <EOL> help = "<STR_LIT>" ) <EOL> p . add_option ( "<STR_LIT>" , default = False , action = "<STR_LIT:store_true>" , <EOL> help = "<STR_LIT>" ) <EOL> p . add_option ( "<STR_LIT>" , default = False , action = "<STR_LIT:store_true>" , <EOL> help = "<STR_LIT>" ) <EOL> opts , args = p . parse_args ( args ) <EOL> if len ( args ) != <NUM_LIT:2> : <EOL> sys . exit ( not p . print_help ( ) ) <EOL> afasta , bfasta = args <EOL> afastan = len ( Fasta ( afasta ) ) <EOL> bfastan = len ( Fasta ( bfasta ) ) <EOL> if afastan == bfastan : <EOL> print green ( "<STR_LIT>" . format ( afastan , bfastan ) ) <EOL> else : <EOL> print red ( "<STR_LIT>" . format ( afastan , bfastan ) ) <EOL> ah = SeqIO . parse ( afasta , "<STR_LIT>" ) <EOL> bh = SeqIO . parse ( bfasta , "<STR_LIT>" ) <EOL> problem_ids = [ ] <EOL> for arec , brec in zip ( ah , bh ) : <EOL> if opts . ignore_stop : <EOL> arec . seq = arec . seq . rstrip ( "<STR_LIT:*>" ) <EOL> brec . seq = brec . seq . rstrip ( "<STR_LIT:*>" ) <EOL> asize , bsize = len ( arec ) , len ( brec ) <EOL> if not opts . quiet : <EOL> print banner ( str ( arec ) , [ str ( brec ) ] ) <EOL> if asize == bsize : <EOL> print green ( "<STR_LIT>" % asize ) <EOL> else : <EOL> print red ( "<STR_LIT>" % ( asize , bsize ) ) <EOL> fd = print_first_difference ( arec , brec , ignore_case = opts . ignore_case , <EOL> ignore_N = opts . ignore_N , rc = opts . rc , report_match = not opts . quiet ) <EOL> if not fd : <EOL> logging . error ( "<STR_LIT>" . format ( arec . id ) ) <EOL> problem_ids . append ( "<STR_LIT:\t>" . join ( str ( x ) for x in ( arec . id , asize , bsize , <EOL> abs ( asize - bsize ) ) ) ) <EOL> if problem_ids : <EOL> print red ( "<STR_LIT>" . format ( len ( problem_ids ) ) ) <EOL> fw = must_open ( "<STR_LIT>" , "<STR_LIT:w>" ) <EOL> print >> fw , "<STR_LIT:\n>" . join ( problem_ids ) <EOL> def hash_fasta ( seq , ignore_case = False , ignore_N = False , ignore_stop = False , checksum = "<STR_LIT>" ) : <EOL> """<STR_LIT>""" <EOL> if ignore_stop : <EOL> seq = seq . rstrip ( "<STR_LIT:*>" ) <EOL> if ignore_case : <EOL> seq = seq . upper ( ) <EOL> if ignore_N : <EOL> if not all ( c . upper ( ) in '<STR_LIT>' for c in seq ) : <EOL> seq = re . sub ( '<STR_LIT:X>' , '<STR_LIT>' , seq ) <EOL> else : <EOL> seq = re . sub ( '<STR_LIT:N>' , '<STR_LIT>' , seq ) <EOL> if checksum == "<STR_LIT>" : <EOL> hashed = md5 ( seq ) . hexdigest ( ) <EOL> elif checksum == "<STR_LIT>" : <EOL> hashed = seguid ( seq ) <EOL> return hashed <EOL> def identical ( args ) : <EOL> """<STR_LIT>""" <EOL> from jcvi . utils . cbook import AutoVivification <EOL> allowed_checksum = [ "<STR_LIT>" , "<STR_LIT>" ] <EOL> p = OptionParser ( identical . __doc__ ) <EOL> p . add_option ( "<STR_LIT>" , default = False , action = "<STR_LIT:store_true>" , <EOL> help = "<STR_LIT>" ) <EOL> p . add_option ( "<STR_LIT>" , default = False , action = "<STR_LIT:store_true>" , <EOL> help = "<STR_LIT>" ) <EOL> p . add_option ( "<STR_LIT>" , default = False , action = "<STR_LIT:store_true>" , <EOL> help = "<STR_LIT>" ) <EOL> p . add_option ( "<STR_LIT>" , default = False , action = "<STR_LIT:store_true>" , <EOL> help = "<STR_LIT>" + "<STR_LIT>" ) <EOL> p . add_option ( "<STR_LIT>" , default = "<STR_LIT>" , choices = allowed_checksum , <EOL> help = "<STR_LIT>" ) <EOL> p . set_outfile ( ) <EOL> opts , args = p . parse_args ( args ) <EOL> if len ( args ) == <NUM_LIT:0> : <EOL> sys . exit ( not p . print_help ( ) ) <EOL> d = AutoVivification ( ) <EOL> files = [ ] <EOL> for fastafile in args : <EOL> f = Fasta ( fastafile ) <EOL> pf = fastafile . rsplit ( "<STR_LIT:.>" , <NUM_LIT:1> ) [ <NUM_LIT:0> ] <EOL> files . append ( pf ) <EOL> logging . debug ( "<STR_LIT>" . format ( fastafile ) ) <EOL> for name , rec in f . iteritems_ordered ( ) : <EOL> seq = re . sub ( '<STR_LIT:U+0020>' , '<STR_LIT>' , str ( rec . seq ) ) <EOL> hashed = hash_fasta ( seq , ignore_case = opts . ignore_case , ignore_N = opts . ignore_N , ignore_stop = opts . ignore_stop , checksum = opts . checksum ) <EOL> if not d [ hashed ] : <EOL> d [ hashed ] [ '<STR_LIT>' ] = seq <EOL> d [ hashed ] [ '<STR_LIT:count>' ] = <NUM_LIT:0> <EOL> if not d [ hashed ] [ '<STR_LIT>' ] [ pf ] : <EOL> d [ hashed ] [ '<STR_LIT>' ] [ pf ] = set ( ) <EOL> d [ hashed ] [ '<STR_LIT>' ] [ pf ] . add ( name ) <EOL> fw = must_open ( opts . outfile , "<STR_LIT:w>" ) <EOL> if opts . output_uniq : <EOL> uniqfile = "<STR_LIT:_>" . join ( files ) + "<STR_LIT>" <EOL> uniqfw = must_open ( uniqfile , "<STR_LIT:w>" ) <EOL> header = "<STR_LIT:\t>" . join ( str ( x ) for x in ( args ) ) <EOL> print >> fw , "<STR_LIT:\t>" . join ( str ( x ) for x in ( "<STR_LIT>" , header ) ) <EOL> for idx , hashed in enumerate ( d . keys ( ) ) : <EOL> line = [ ] <EOL> line . append ( "<STR_LIT>" . format ( idx ) ) <EOL> for fastafile in files : <EOL> if fastafile in d [ hashed ] [ '<STR_LIT>' ] . keys ( ) : <EOL> line . append ( "<STR_LIT:U+002C>" . join ( d [ hashed ] [ '<STR_LIT>' ] [ fastafile ] ) ) <EOL> if opts . output_uniq : <EOL> d [ hashed ] [ '<STR_LIT:count>' ] += len ( d [ hashed ] [ '<STR_LIT>' ] [ fastafile ] ) <EOL> else : <EOL> line . append ( "<STR_LIT>" ) <EOL> print >> fw , "<STR_LIT:\t>" . join ( line ) <EOL> if opts . output_uniq : <EOL> seqid = "<STR_LIT:\t>" . join ( str ( x ) for x in ( "<STR_LIT>" . format ( idx ) , d [ hashed ] [ '<STR_LIT:count>' ] ) ) <EOL> rec = SeqRecord ( Seq ( d [ hashed ] [ '<STR_LIT>' ] ) , id = seqid , description = "<STR_LIT>" ) <EOL> SeqIO . write ( [ rec ] , uniqfw , "<STR_LIT>" ) <EOL> fw . close ( ) <EOL> if opts . output_uniq : <EOL> logging . debug ( "<STR_LIT>" . format ( uniqfile ) ) <EOL> uniqfw . close ( ) <EOL> QUALSUFFIX = "<STR_LIT>" <EOL> def get_qual ( fastafile , suffix = QUALSUFFIX , check = True ) : <EOL> """<STR_LIT>""" <EOL> qualfile1 = fastafile . rsplit ( "<STR_LIT:.>" , <NUM_LIT:1> ) [ <NUM_LIT:0> ] + suffix <EOL> qualfile2 = fastafile + suffix <EOL> if check : <EOL> if op . exists ( qualfile1 ) : <EOL> logging . debug ( "<STR_LIT>" . format ( qualfile1 ) ) <EOL> return qualfile1 <EOL> elif op . exists ( qualfile2 ) : <EOL> logging . debug ( "<STR_LIT>" . format ( qualfile2 ) ) <EOL> return qualfile2 <EOL> else : <EOL> logging . warning ( "<STR_LIT>" ) <EOL> return None <EOL> return qualfile1 <EOL> def some ( args ) : <EOL> """<STR_LIT>""" <EOL> p = OptionParser ( some . __doc__ ) <EOL> p . add_option ( "<STR_LIT>" , default = False , action = "<STR_LIT:store_true>" , <EOL> help = "<STR_LIT>" ) <EOL> p . add_option ( "<STR_LIT>" , default = False , action = "<STR_LIT:store_true>" , <EOL> help = "<STR_LIT>" ) <EOL> opts , args = p . parse_args ( args ) <EOL> if len ( args ) != <NUM_LIT:3> : <EOL> sys . exit ( p . print_help ( ) ) <EOL> fastafile , listfile , outfastafile = args <EOL> outfastahandle = must_open ( outfastafile , "<STR_LIT:w>" ) <EOL> qualfile = get_qual ( fastafile ) <EOL> names = set ( x . strip ( ) for x in open ( listfile ) ) <EOL> if qualfile : <EOL> outqualfile = outfastafile + "<STR_LIT>" <EOL> outqualhandle = open ( outqualfile , "<STR_LIT:w>" ) <EOL> parser = iter_fasta_qual ( fastafile , qualfile ) <EOL> else : <EOL> parser = SeqIO . parse ( fastafile , "<STR_LIT>" ) <EOL> num_records = <NUM_LIT:0> <EOL> for rec in parser : <EOL> name = rec . id <EOL> if opts . uniprot : <EOL> name = name . split ( "<STR_LIT:|>" ) [ - <NUM_LIT:1> ] <EOL> if opts . exclude : <EOL> if name in names : <EOL> continue <EOL> else : <EOL> if name not in names : <EOL> continue <EOL> SeqIO . write ( [ rec ] , outfastahandle , "<STR_LIT>" ) <EOL> if qualfile : <EOL> SeqIO . write ( [ rec ] , outqualhandle , "<STR_LIT>" ) <EOL> num_records += <NUM_LIT:1> <EOL> logging . debug ( "<STR_LIT>" % ( num_records , outfastafile ) ) <EOL> def fastq ( args ) : <EOL> """<STR_LIT>""" <EOL> from jcvi . formats . fastq import FastqLite <EOL> p = OptionParser ( fastq . __doc__ ) <EOL> p . add_option ( "<STR_LIT>" , type = "<STR_LIT:int>" , <EOL> help = "<STR_LIT>" ) <EOL> opts , args = p . parse_args ( args ) <EOL> if len ( args ) != <NUM_LIT:1> : <EOL> sys . exit ( not p . print_help ( ) ) <EOL> fastafile , = args <EOL> fastqfile = fastafile . rsplit ( "<STR_LIT:.>" , <NUM_LIT:1> ) [ <NUM_LIT:0> ] + "<STR_LIT>" <EOL> fastqhandle = open ( fastqfile , "<STR_LIT:w>" ) <EOL> num_records = <NUM_LIT:0> <EOL> if opts . qv is not None : <EOL> qv = chr ( ord ( '<STR_LIT:!>' ) + opts . qv ) <EOL> logging . debug ( "<STR_LIT>" . format ( qv , opts . qv ) ) <EOL> else : <EOL> qv = None <EOL> if qv : <EOL> f = Fasta ( fastafile , lazy = True ) <EOL> for name , rec in f . iteritems_ordered ( ) : <EOL> r = FastqLite ( "<STR_LIT:@>" + name , str ( rec . seq ) . upper ( ) , qv * len ( rec . seq ) ) <EOL> print >> fastqhandle , r <EOL> num_records += <NUM_LIT:1> <EOL> else : <EOL> qualfile = get_qual ( fastafile ) <EOL> for rec in iter_fasta_qual ( fastafile , qualfile ) : <EOL> SeqIO . write ( [ rec ] , fastqhandle , "<STR_LIT>" ) <EOL> num_records += <NUM_LIT:1> <EOL> fastqhandle . close ( ) <EOL> logging . debug ( "<STR_LIT>" % ( num_records , fastqfile ) ) <EOL> def pair ( args ) : <EOL> """<STR_LIT>""" <EOL> p = OptionParser ( pair . __doc__ ) <EOL> p . set_sep ( sep = None , help = "<STR_LIT>" + "<STR_LIT>" ) <EOL> p . add_option ( "<STR_LIT>" , dest = "<STR_LIT>" , default = False , action = "<STR_LIT:store_true>" , <EOL> help = "<STR_LIT>" ) <EOL> opts , args = p . parse_args ( args ) <EOL> if len ( args ) != <NUM_LIT:1> : <EOL> sys . exit ( p . print_help ( ) ) <EOL> fastafile , = args <EOL> qualfile = get_qual ( fastafile ) <EOL> prefix = fastafile . rsplit ( "<STR_LIT:.>" , <NUM_LIT:1> ) [ <NUM_LIT:0> ] <EOL> pairsfile = prefix + "<STR_LIT>" <EOL> fragsfile = prefix + "<STR_LIT>" <EOL> pairsfw = open ( pairsfile , "<STR_LIT:w>" ) <EOL> fragsfw = open ( fragsfile , "<STR_LIT:w>" ) <EOL> if opts . matepairs : <EOL> matepairsfile = prefix + "<STR_LIT>" <EOL> matepairsfw = open ( matepairsfile , "<STR_LIT:w>" ) <EOL> if qualfile : <EOL> pairsqualfile = pairsfile + "<STR_LIT>" <EOL> pairsqualhandle = open ( pairsqualfile , "<STR_LIT:w>" ) <EOL> fragsqualfile = fragsfile + "<STR_LIT>" <EOL> fragsqualhandle = open ( fragsqualfile , "<STR_LIT:w>" ) <EOL> f = Fasta ( fastafile ) <EOL> if qualfile : <EOL> q = SeqIO . index ( qualfile , "<STR_LIT>" ) <EOL> all_keys = list ( f . iterkeys ( ) ) <EOL> all_keys . sort ( ) <EOL> sep = opts . sep <EOL> if sep : <EOL> key_fun = lambda x : x . split ( sep , <NUM_LIT:1> ) [ <NUM_LIT:0> ] <EOL> else : <EOL> key_fun = lambda x : x [ : - <NUM_LIT:1> ] <EOL> for key , variants in groupby ( all_keys , key = key_fun ) : <EOL> variants = list ( variants ) <EOL> paired = ( len ( variants ) == <NUM_LIT:2> ) <EOL> if paired and opts . matepairs : <EOL> print >> matepairsfw , "<STR_LIT:\t>" . join ( ( "<STR_LIT>" % key , "<STR_LIT>" % key ) ) <EOL> fw = pairsfw if paired else fragsfw <EOL> if qualfile : <EOL> qualfw = pairsqualhandle if paired else fragsqualhandle <EOL> for i , var in enumerate ( variants ) : <EOL> rec = f [ var ] <EOL> if qualfile : <EOL> recqual = q [ var ] <EOL> newid = "<STR_LIT>" % ( key , i + <NUM_LIT:1> ) <EOL> rec . id = newid <EOL> rec . description = "<STR_LIT>" <EOL> SeqIO . write ( [ rec ] , fw , "<STR_LIT>" ) <EOL> if qualfile : <EOL> recqual . id = newid <EOL> recqual . description = "<STR_LIT>" <EOL> SeqIO . write ( [ recqual ] , qualfw , "<STR_LIT>" ) <EOL> logging . debug ( "<STR_LIT>" % ( pairsfile , fragsfile ) ) <EOL> if opts . matepairs : <EOL> logging . debug ( "<STR_LIT>" % matepairsfile ) <EOL> def pairinplace ( args ) : <EOL> """<STR_LIT>""" <EOL> from jcvi . utils . iter import pairwise <EOL> p = OptionParser ( pairinplace . __doc__ ) <EOL> p . add_option ( "<STR_LIT>" , dest = "<STR_LIT>" , default = <NUM_LIT:1> , type = "<STR_LIT:int>" , <EOL> help = "<STR_LIT>" ) <EOL> opts , args = p . parse_args ( args ) <EOL> if len ( args ) != <NUM_LIT:1> : <EOL> sys . exit ( not p . print_help ( ) ) <EOL> fastafile , = args <EOL> base = op . basename ( fastafile ) . split ( "<STR_LIT:.>" ) [ <NUM_LIT:0> ] <EOL> frags = base + "<STR_LIT>" <EOL> pairs = base + "<STR_LIT>" <EOL> if fastafile . endswith ( "<STR_LIT>" ) : <EOL> frags += "<STR_LIT>" <EOL> pairs += "<STR_LIT>" <EOL> fragsfw = must_open ( frags , "<STR_LIT:w>" ) <EOL> pairsfw = must_open ( pairs , "<STR_LIT:w>" ) <EOL> N = opts . rclip <EOL> strip_name = lambda x : x [ : - N ] if N else str <EOL> skipflag = False <EOL> fastaiter = SeqIO . parse ( fastafile , "<STR_LIT>" ) <EOL> for a , b in pairwise ( fastaiter ) : <EOL> aid , bid = [ strip_name ( x ) for x in ( a . id , b . id ) ] <EOL> if skipflag : <EOL> skipflag = False <EOL> continue <EOL> if aid == bid : <EOL> SeqIO . write ( [ a , b ] , pairsfw , "<STR_LIT>" ) <EOL> skipflag = True <EOL> else : <EOL> SeqIO . write ( [ a ] , fragsfw , "<STR_LIT>" ) <EOL> if not skipflag : <EOL> SeqIO . write ( [ a ] , fragsfw , "<STR_LIT>" ) <EOL> logging . debug ( "<STR_LIT>" % ( pairs , frags ) ) <EOL> def extract ( args ) : <EOL> """<STR_LIT>""" <EOL> p = OptionParser ( extract . __doc__ ) <EOL> p . add_option ( '<STR_LIT>' , default = False , action = "<STR_LIT:store_true>" , <EOL> help = "<STR_LIT>" ) <EOL> p . add_option ( '<STR_LIT>' , default = False , action = "<STR_LIT:store_true>" , <EOL> help = "<STR_LIT>" ) <EOL> p . add_option ( '<STR_LIT>' , default = False , action = "<STR_LIT:store_true>" , <EOL> help = "<STR_LIT>" ) <EOL> p . add_option ( '<STR_LIT>' , default = None , <EOL> help = "<STR_LIT>" ) <EOL> p . set_outfile ( ) <EOL> opts , args = p . parse_args ( args ) <EOL> if len ( args ) == <NUM_LIT:2> : <EOL> fastafile , query = args <EOL> elif len ( args ) == <NUM_LIT:1> and opts . bed : <EOL> fastafile , = args <EOL> bedaccns = Bed ( opts . bed ) . accns <EOL> else : <EOL> sys . exit ( p . print_help ( ) ) <EOL> if opts . bed : <EOL> fw = must_open ( opts . outfile , "<STR_LIT:w>" ) <EOL> f = Fasta ( fastafile ) <EOL> for accn in bedaccns : <EOL> try : <EOL> rec = f [ accn ] <EOL> except : <EOL> logging . error ( "<STR_LIT>" . format ( accn , fastafile ) ) <EOL> continue <EOL> SeqIO . write ( [ rec ] , fw , "<STR_LIT>" ) <EOL> return fw . name <EOL> atoms = query . split ( "<STR_LIT::>" ) <EOL> key = atoms [ <NUM_LIT:0> ] <EOL> assert len ( atoms ) <= <NUM_LIT:3> , "<STR_LIT>" <EOL> pos = "<STR_LIT>" <EOL> if len ( atoms ) in ( <NUM_LIT:2> , <NUM_LIT:3> ) : <EOL> pos = atoms [ <NUM_LIT:1> ] <EOL> strand = "<STR_LIT:+>" <EOL> if len ( atoms ) == <NUM_LIT:3> : <EOL> strand = atoms [ <NUM_LIT:2> ] <EOL> assert strand in ( '<STR_LIT:+>' , '<STR_LIT:->' ) , "<STR_LIT>" <EOL> feature = dict ( chr = key ) <EOL> if "<STR_LIT:->" in pos : <EOL> start , stop = pos . split ( "<STR_LIT:->" ) <EOL> try : <EOL> start , stop = int ( start ) , int ( stop ) <EOL> except ValueError as e : <EOL> logging . error ( e ) <EOL> sys . exit ( p . print_help ( ) ) <EOL> feature [ "<STR_LIT:start>" ] = start <EOL> feature [ "<STR_LIT>" ] = stop <EOL> else : <EOL> start , stop = None , None <EOL> assert start < stop or None in ( start , stop ) , "<STR_LIT>" . format ( start , stop ) <EOL> feature [ "<STR_LIT>" ] = strand <EOL> include , exclude = opts . include , opts . exclude <EOL> assert not ( include and exclude ) , "<STR_LIT>" "<STR_LIT>" <EOL> fw = must_open ( opts . outfile , "<STR_LIT:w>" ) <EOL> if include or exclude : <EOL> f = Fasta ( fastafile , lazy = True ) <EOL> fi = f . iteritems_ordered if opts . idonly else f . iterdescriptions_ordered <EOL> for k , rec in fi ( ) : <EOL> if include and key not in k : <EOL> continue <EOL> if exclude and key in k : <EOL> continue <EOL> seq = Fasta . subseq ( rec , start , stop , strand ) <EOL> newid = rec . id <EOL> if start is not None : <EOL> newid += "<STR_LIT>" . format ( start , stop , strand ) <EOL> rec = SeqRecord ( seq , id = newid , description = k ) <EOL> SeqIO . write ( [ rec ] , fw , "<STR_LIT>" ) <EOL> else : <EOL> f = Fasta ( fastafile ) <EOL> try : <EOL> seq = f . sequence ( feature , asstring = False ) <EOL> except AssertionError as e : <EOL> logging . error ( e ) <EOL> return <EOL> rec = SeqRecord ( seq , id = query , description = "<STR_LIT>" ) <EOL> SeqIO . write ( [ rec ] , fw , "<STR_LIT>" ) <EOL> return fw . name <EOL> def _uniq_rec ( fastafile , seq = False ) : <EOL> """<STR_LIT>""" <EOL> seen = set ( ) <EOL> for rec in SeqIO . parse ( fastafile , "<STR_LIT>" ) : <EOL> name = str ( rec . seq ) if seq else rec . id <EOL> if name in seen : <EOL> logging . debug ( "<STR_LIT>" . format ( rec . id ) ) <EOL> continue <EOL> seen . add ( name ) <EOL> yield rec <EOL> def uniq ( args ) : <EOL> """<STR_LIT>""" <EOL> p = OptionParser ( uniq . __doc__ ) <EOL> p . add_option ( "<STR_LIT>" , default = False , action = "<STR_LIT:store_true>" , <EOL> help = "<STR_LIT>" ) <EOL> p . add_option ( "<STR_LIT>" , "<STR_LIT>" , dest = "<STR_LIT>" , <EOL> action = "<STR_LIT:store_true>" , default = False , <EOL> help = "<STR_LIT>" ) <EOL> opts , args = p . parse_args ( args ) <EOL> if len ( args ) != <NUM_LIT:2> : <EOL> sys . exit ( p . print_help ( ) ) <EOL> fastafile , uniqfastafile = args <EOL> fw = must_open ( uniqfastafile , "<STR_LIT:w>" ) <EOL> seq = opts . seq <EOL> for rec in _uniq_rec ( fastafile , seq = seq ) : <EOL> if opts . trimname : <EOL> rec . description = "<STR_LIT>" <EOL> SeqIO . write ( [ rec ] , fw , "<STR_LIT>" ) <EOL> def random ( args ) : <EOL> """<STR_LIT>""" <EOL> from random import sample <EOL> p = OptionParser ( random . __doc__ ) <EOL> opts , args = p . parse_args ( args ) <EOL> if len ( args ) != <NUM_LIT:2> : <EOL> sys . exit ( not p . print_help ( ) ) <EOL> fastafile , N = args <EOL> N = int ( N ) <EOL> assert N > <NUM_LIT:0> <EOL> f = Fasta ( fastafile ) <EOL> fw = must_open ( "<STR_LIT>" , "<STR_LIT:w>" ) <EOL> for key in sample ( f . keys ( ) , N ) : <EOL> rec = f [ key ] <EOL> SeqIO . write ( [ rec ] , fw , "<STR_LIT>" ) <EOL> fw . close ( ) <EOL> XQUAL = - <NUM_LIT:1000> <EOL> NQUAL = <NUM_LIT:5> <EOL> QUAL = <NUM_LIT:10> <EOL> OKQUAL = <NUM_LIT:15> <EOL> def modify_qual ( rec ) : <EOL> qv = rec . letter_annotations [ '<STR_LIT>' ] <EOL> for i , ( s , q ) in enumerate ( zip ( rec . seq , qv ) ) : <EOL> if s == '<STR_LIT:X>' or s == '<STR_LIT:x>' : <EOL> qv [ i ] = XQUAL <EOL> if s == '<STR_LIT:N>' or s == '<STR_LIT:x>' : <EOL> qv [ i ] = NQUAL <EOL> return rec <EOL> def make_qual ( fastafile , score = OKQUAL ) : <EOL> logging . warning ( "<STR_LIT>" . format ( score ) ) <EOL> qualfile = fastafile . rsplit ( "<STR_LIT:.>" , <NUM_LIT:1> ) [ <NUM_LIT:0> ] + "<STR_LIT>" <EOL> fw = open ( qualfile , "<STR_LIT:w>" ) <EOL> fasta = Fasta ( fastafile , lazy = True ) <EOL> score = str ( score ) + "<STR_LIT:U+0020>" <EOL> for entry , size in fasta . itersizes_ordered ( ) : <EOL> print >> fw , "<STR_LIT:>>" + entry <EOL> print >> fw , score * size <EOL> fw . close ( ) <EOL> return qualfile <EOL> def iter_fasta_qual ( fastafile , qualfile , defaultqual = OKQUAL , modify = False ) : <EOL> """<STR_LIT>""" <EOL> from Bio . SeqIO . QualityIO import PairedFastaQualIterator <EOL> if not qualfile : <EOL> qualfile = make_qual ( fastafile , score = defaultqual ) <EOL> rec_iter = PairedFastaQualIterator ( open ( fastafile ) , open ( qualfile ) ) <EOL> for rec in rec_iter : <EOL> yield rec if not modify else modify_qual ( rec ) <EOL> def write_fasta_qual ( rec , fastahandle , qualhandle ) : <EOL> if fastahandle : <EOL> SeqIO . write ( [ rec ] , fastahandle , "<STR_LIT>" ) <EOL> if qualhandle : <EOL> SeqIO . write ( [ rec ] , qualhandle , "<STR_LIT>" ) <EOL> def trim ( args ) : <EOL> """<STR_LIT>""" <EOL> from jcvi . algorithms . maxsum import max_sum <EOL> p = OptionParser ( trim . __doc__ ) <EOL> p . add_option ( "<STR_LIT:-c>" , dest = "<STR_LIT>" , type = "<STR_LIT:int>" , default = <NUM_LIT:64> , <EOL> help = "<STR_LIT>" ) <EOL> p . add_option ( "<STR_LIT>" , dest = "<STR_LIT>" , default = QUAL , <EOL> help = "<STR_LIT>" ) <EOL> opts , args = p . parse_args ( args ) <EOL> if len ( args ) != <NUM_LIT:2> : <EOL> sys . exit ( p . print_help ( ) ) <EOL> fastafile , newfastafile = args <EOL> qualfile = get_qual ( fastafile ) <EOL> newqualfile = get_qual ( newfastafile , check = False ) <EOL> logging . debug ( "<STR_LIT>" % ( fastafile , newfastafile ) ) <EOL> fw = must_open ( newfastafile , "<STR_LIT:w>" ) <EOL> fw_qual = open ( newqualfile , "<STR_LIT:w>" ) <EOL> dropped = trimmed = <NUM_LIT:0> <EOL> for rec in iter_fasta_qual ( fastafile , qualfile , modify = True ) : <EOL> qv = [ x - opts . score for x in rec . letter_annotations [ "<STR_LIT>" ] ] <EOL> msum , trim_start , trim_end = max_sum ( qv ) <EOL> score = trim_end - trim_start + <NUM_LIT:1> <EOL> if score < opts . min_length : <EOL> dropped += <NUM_LIT:1> <EOL> continue <EOL> if score < len ( rec ) : <EOL> trimmed += <NUM_LIT:1> <EOL> rec = rec [ trim_start : trim_end + <NUM_LIT:1> ] <EOL> write_fasta_qual ( rec , fw , fw_qual ) <EOL> print >> sys . stderr , "<STR_LIT>" % trimmed <EOL> print >> sys . stderr , "<STR_LIT>" % ( dropped , opts . min_length ) <EOL> fw . close ( ) <EOL> fw_qual . close ( ) <EOL> def sequin ( args ) : <EOL> """<STR_LIT>""" <EOL> p = OptionParser ( sequin . __doc__ ) <EOL> p . add_option ( "<STR_LIT>" , default = <NUM_LIT:100> , type = "<STR_LIT:int>" , <EOL> help = "<STR_LIT>" ) <EOL> p . add_option ( "<STR_LIT>" , default = None , <EOL> help = "<STR_LIT>" ) <EOL> p . add_option ( "<STR_LIT>" , default = None , <EOL> help = "<STR_LIT>" ) <EOL> p . add_option ( "<STR_LIT>" , default = None , <EOL> help = "<STR_LIT>" ) <EOL> p . set_mingap ( default = <NUM_LIT:100> ) <EOL> opts , args = p . parse_args ( args ) <EOL> if len ( args ) != <NUM_LIT:1> : <EOL> sys . exit ( not p . print_help ( ) ) <EOL> inputfasta , = args <EOL> unk = opts . unk <EOL> outputfasta = inputfasta . rsplit ( "<STR_LIT:.>" , <NUM_LIT:1> ) [ <NUM_LIT:0> ] + "<STR_LIT>" <EOL> rec = SeqIO . parse ( must_open ( inputfasta ) , "<STR_LIT>" ) . next ( ) <EOL> seq = "<STR_LIT>" <EOL> unknowns , knowns = <NUM_LIT:0> , <NUM_LIT:0> <EOL> for gap , gap_group in groupby ( rec . seq , lambda x : x . upper ( ) == '<STR_LIT:N>' ) : <EOL> subseq = "<STR_LIT>" . join ( gap_group ) <EOL> if gap : <EOL> gap_length = len ( subseq ) <EOL> if gap_length == unk : <EOL> subseq = "<STR_LIT>" . format ( unk ) <EOL> unknowns += <NUM_LIT:1> <EOL> elif gap_length >= opts . mingap : <EOL> subseq = "<STR_LIT>" . format ( gap_length ) <EOL> knowns += <NUM_LIT:1> <EOL> seq += subseq <EOL> fw = must_open ( outputfasta , "<STR_LIT:w>" ) <EOL> id = opts . newid or rec . id <EOL> fastaheader = "<STR_LIT>" . format ( id ) <EOL> if opts . chromosome : <EOL> fastaheader += "<STR_LIT>" . format ( opts . chromosome ) <EOL> if opts . clone : <EOL> fastaheader += "<STR_LIT>" . format ( opts . clone ) <EOL> print >> fw , fastaheader <EOL> print >> fw , seq <EOL> fw . close ( ) <EOL> logging . debug ( "<STR_LIT>" . format ( outputfasta , unknowns , knowns ) ) <EOL> return outputfasta , unknowns + knowns <EOL> def remove_small_components ( rec , minlen ) : <EOL> newseq = [ ] <EOL> removed = <NUM_LIT:0> <EOL> for gap , seq in groupby ( rec . seq , lambda x : x . upper ( ) == '<STR_LIT:N>' ) : <EOL> seq = "<STR_LIT>" . join ( seq ) <EOL> seqlen = len ( seq ) <EOL> if not gap and seqlen < minlen : <EOL> seq = seqlen * "<STR_LIT:N>" <EOL> logging . debug ( "<STR_LIT>" . format ( seqlen , rec . name ) ) <EOL> removed += seqlen <EOL> newseq . append ( seq ) <EOL> rec . seq = Seq ( "<STR_LIT>" . join ( newseq ) ) <EOL> return removed <EOL> def trim_terminal_Ns ( rec ) : <EOL> rec . seq = rec . seq . strip ( "<STR_LIT:N>" ) <EOL> def normalize_gaps ( rec , gapsize ) : <EOL> newseq = [ ] <EOL> normalized = <NUM_LIT:0> <EOL> NN = gapsize * "<STR_LIT:N>" <EOL> for gap , seq in groupby ( rec . seq , lambda x : x . upper ( ) == '<STR_LIT:N>' ) : <EOL> seq = "<STR_LIT>" . join ( seq ) <EOL> if gap : <EOL> seq = NN <EOL> normalized += <NUM_LIT:1> <EOL> newseq . append ( seq ) <EOL> rec . seq = Seq ( "<STR_LIT>" . join ( newseq ) ) <EOL> return normalized <EOL> def tidy ( args ) : <EOL> """<STR_LIT>""" <EOL> p = OptionParser ( tidy . __doc__ ) <EOL> p . add_option ( "<STR_LIT>" , dest = "<STR_LIT>" , default = <NUM_LIT:0> , type = "<STR_LIT:int>" , <EOL> help = "<STR_LIT>" ) <EOL> p . add_option ( "<STR_LIT>" , dest = "<STR_LIT>" , default = <NUM_LIT:100> , type = "<STR_LIT:int>" , <EOL> help = "<STR_LIT>" ) <EOL> opts , args = p . parse_args ( args ) <EOL> if len ( args ) != <NUM_LIT:1> : <EOL> sys . exit ( not p . print_help ( ) ) <EOL> fastafile , = args <EOL> gapsize = opts . gapsize <EOL> minlen = opts . minlen <EOL> tidyfastafile = fastafile . rsplit ( "<STR_LIT:.>" , <NUM_LIT:1> ) [ <NUM_LIT:0> ] + "<STR_LIT>" <EOL> fw = must_open ( tidyfastafile , "<STR_LIT:w>" ) <EOL> removed = normalized = <NUM_LIT:0> <EOL> fasta = Fasta ( fastafile , lazy = True ) <EOL> for name , rec in fasta . iteritems_ordered ( ) : <EOL> rec . seq = rec . seq . upper ( ) <EOL> if minlen : <EOL> removed += remove_small_components ( rec , minlen ) <EOL> trim_terminal_Ns ( rec ) <EOL> if gapsize : <EOL> normalized += normalize_gaps ( rec , gapsize ) <EOL> if len ( rec ) == <NUM_LIT:0> : <EOL> logging . debug ( "<STR_LIT>" . format ( rec . id ) ) <EOL> continue <EOL> SeqIO . write ( [ rec ] , fw , "<STR_LIT>" ) <EOL> if removed : <EOL> logging . debug ( "<STR_LIT>" . format ( removed ) ) <EOL> if normalized : <EOL> logging . debug ( "<STR_LIT>" . format ( normalized ) ) <EOL> logging . debug ( "<STR_LIT>" . format ( tidyfastafile ) ) <EOL> fw . close ( ) <EOL> return tidyfastafile <EOL> def write_gaps_worker ( rec ) : <EOL> start = <NUM_LIT:0> <EOL> seq = rec . seq . upper ( ) <EOL> output = [ ] <EOL> for gap , seq in groupby ( seq , lambda x : x == '<STR_LIT:N>' ) : <EOL> seq = "<STR_LIT>" . join ( seq ) <EOL> current_length = len ( seq ) <EOL> object_beg = start + <NUM_LIT:1> <EOL> object_end = start + current_length <EOL> if gap : <EOL> s = "<STR_LIT:\t>" . join ( str ( x ) for x in ( rec . id , <EOL> object_beg - <NUM_LIT:1> , object_end ) ) <EOL> output . append ( s ) <EOL> start += current_length <EOL> return "<STR_LIT:\n>" . join ( output ) <EOL> def write_gaps_bed ( inputfasta , prefix , mingap , cpus ) : <EOL> from jcvi . apps . grid import WriteJobs <EOL> from jcvi . formats . bed import sort <EOL> bedfile = prefix + "<STR_LIT>" <EOL> f = Fasta ( inputfasta ) <EOL> recs = list ( rec for k , rec in f . iteritems ( ) ) <EOL> pool = WriteJobs ( write_gaps_worker , recs , bedfile , cpus = cpus ) <EOL> pool . run ( ) <EOL> sort ( [ bedfile , "<STR_LIT>" ] ) <EOL> bed = Bed ( bedfile ) <EOL> nbedfile = prefix + "<STR_LIT>" . format ( mingap ) <EOL> gapnum = <NUM_LIT:0> <EOL> fw = open ( nbedfile , "<STR_LIT:w>" ) <EOL> for b in bed : <EOL> if b . span < mingap : <EOL> continue <EOL> gapnum += <NUM_LIT:1> <EOL> gapname = "<STR_LIT>" . format ( gapnum ) <EOL> print >> fw , "<STR_LIT:\t>" . join ( str ( x ) for x in ( b , gapname , b . span ) ) <EOL> shutil . move ( nbedfile , bedfile ) <EOL> logging . debug ( "<STR_LIT>" . format ( mingap , bedfile ) ) <EOL> def gaps ( args ) : <EOL> """<STR_LIT>""" <EOL> from jcvi . formats . sizes import agp <EOL> from jcvi . formats . agp import mask , build <EOL> p = OptionParser ( gaps . __doc__ ) <EOL> p . add_option ( "<STR_LIT>" , default = False , action = "<STR_LIT:store_true>" , <EOL> help = "<STR_LIT>" ) <EOL> p . set_mingap ( default = <NUM_LIT:100> ) <EOL> p . set_cpus ( ) <EOL> opts , args = p . parse_args ( args ) <EOL> if len ( args ) != <NUM_LIT:1> : <EOL> sys . exit ( not p . print_help ( ) ) <EOL> inputfasta , = args <EOL> mingap = opts . mingap <EOL> split = opts . split <EOL> prefix = inputfasta . rsplit ( "<STR_LIT:.>" , <NUM_LIT:1> ) [ <NUM_LIT:0> ] <EOL> bedfile = prefix + "<STR_LIT>" <EOL> if need_update ( inputfasta , bedfile ) : <EOL> write_gaps_bed ( inputfasta , prefix , mingap , opts . cpus ) <EOL> if split : <EOL> splitfile = prefix + "<STR_LIT>" <EOL> oagpfile = prefix + "<STR_LIT>" <EOL> cagpfile = prefix + "<STR_LIT>" <EOL> if need_update ( ( inputfasta , bedfile ) , splitfile ) : <EOL> sizesagpfile = agp ( [ inputfasta ] ) <EOL> maskedagpfile = mask ( [ sizesagpfile , bedfile , "<STR_LIT>" ] ) <EOL> shutil . move ( maskedagpfile , oagpfile ) <EOL> logging . debug ( "<STR_LIT>" . format ( oagpfile ) ) <EOL> maskedagpfile = mask ( [ sizesagpfile , bedfile , "<STR_LIT>" ] ) <EOL> shutil . move ( maskedagpfile , cagpfile ) <EOL> logging . debug ( "<STR_LIT>" . format ( cagpfile ) ) <EOL> build ( [ oagpfile , inputfasta , splitfile ] ) <EOL> os . remove ( sizesagpfile ) <EOL> return splitfile , oagpfile , cagpfile <EOL> if __name__ == '<STR_LIT:__main__>' : <EOL> main ( ) </s>
