creates GFE from HLA sequence and locus

        :param locus: string containing HLA locus.
        :param sequence: string containing sequence data.

        :return: GFEobject.
493875
def _breakup_gfe(self, gfe):
        
        [locus, feature_accessions] = gfe.split("w")
        accessions = feature_accessions.split("-")
        i = 0
        features = {}
        for feature_rank in self.structures[locus]:
            accession = accessions[i]
            features.update({feature_rank: accession})
            i += 1

        return(features)
394493
def sequence_type(seq):
    '''

    '''
    if isinstance(seq, coral.DNA):
        material = 'dna'
    elif isinstance(seq, coral.RNA):
        material = 'rna'
    elif isinstance(seq, coral.Peptide):
        material = 'peptide'
    else:
        raise ValueError('Input was not a recognized coral.sequence object.')
    return material
284315
def rev_comp( seq, molecule='dna' ):
	
	
	if molecule == 'dna':
		nuc_dict = { "A":"T", "B":"V", "C":"G", "D":"H", "G":"C", "H":"D", "K":"M", "M":"K", "N":"N", "R":"Y", "S":"S", "T":"A", "V":"B", "W":"W", "Y":"R" }
	elif molecule == 'rna':
		nuc_dict = { "A":"U", "B":"V", "C":"G", "D":"H", "G":"C", "H":"D", "K":"M", "M":"K", "N":"N", "R":"Y", "S":"S", "U":"A", "V":"B", "W":"W", "Y":"R" }
	else:
		raise ValueError( "rev_comp requires molecule to be dna or rna" )
	
	if not isinstance( seq, six.string_types ):
		raise TypeError( "seq must be a string!" )
	
	return ''.join( [ nuc_dict[c] for c in seq.upper()[::-1] ] )
200390
def from_sequence(chain, list_of_residues, sequence_type = None):
        ''''''
        s = Sequence(sequence_type)
        count = 1
        for ResidueAA in list_of_residues:
            s.add(Residue(chain, count, ResidueAA, sequence_type))
            count += 1
        return s
670565
def create_dexseq_annotation(gff, count_file):
    
    out_file = count_file + ".ann"
    if file_exists(out_file):
        return out_file
    with file_transaction(out_file) as tx_out:
        with open(tx_out, 'w') as out_handle:
            with open(gff) as in_handle:
                for line in in_handle:
                    cols = line.strip().split("\t")
                    if cols[2] == "exonic_part":
                        exon = [f for f in cols[8].split(";") if f.strip().startswith("exonic_part_number")]
                        gene = [f for f in cols[8].split(";") if f.strip().startswith("gene_id")]
                        exon = exon[0].replace("\"", "").split()[1]
                        gene = gene[0].replace("\"", "").split()[1]
                        length = int(cols[4]) - int(cols[3]) + 1
                        line = "%s:%s\t%s\t%s\t%s\t%s\t%s\t%s\n" % (gene, exon, gene,
                                                                    cols[0], cols[3],
                                                                    cols[4],
                                                                    length, cols[6])
                        out_handle.write(line)
789791
def build(cls, seqs: Iterable[int], uid: bool = False) -> 'SequenceSet':
        
        seqs_list = sorted(set(seqs))
        groups: List[Union[int, Tuple[int, int]]] = []
        group: Union[int, Tuple[int, int]] = seqs_list[0]
        for i in range(1, len(seqs_list)):
            group_i = seqs_list[i]
            if isinstance(group, int):
                if group_i == group + 1:
                    group = (group, group_i)
                else:
                    groups.append(group)
                    group = group_i
            elif isinstance(group, tuple):
                if group_i == group[1] + 1:
                    group = (group[0], group_i)
                else:
                    groups.append(group)
                    group = group_i
        groups.append(group)
        return SequenceSet(groups, uid)
1057104
def create_seq(self, ):
        
        name = self.name_le.text()
        desc = self.desc_pte.toPlainText()
        try:
            seq = djadapter.models.Sequence(name=name, project=self._project, description=desc)
            seq.save()
            self.sequence = seq
            self.accept()
        except:
            log.exception("Could not create new sequence")
629258
def copy(self):
        '''

        '''
        return type(self)(self.name, self.start, self.stop, self.feature_type,
                          gene=self.gene, locus_tag=self.locus_tag,
                          qualifiers=self.qualifiers, strand=self.strand)
183462
def from_list(cls, gene_ontology, l):
         
        assert isinstance(gene_ontology, GeneOntology)
        assert isinstance(l, list)

        assert len(l) == 17

        go_term = gene_ontology[l[4]]

        qualifier = l[3] or []
        with_from = l[7] or None
        db_name = l[9] or None
        db_syn = l[10] or []
        ext = l[15] or []
        product_id = l[16] or None

        annotation = cls(
            db=l[0],
            db_id=l[1],
            db_symbol=l[2],
            go_term=go_term,
            db_ref=l[5],
            ev_code=l[6],
            db_type=l[11],
            taxon=l[12],
            date=l[13],
            assigned_by=l[14],

            qualifier=qualifier,
            with_from=with_from,
            db_name=db_name,
            db_syn=db_syn,
            ext=ext,
            product_id=product_id,
        )
        return annotation
1029772
def sequencetyper(self):
        
        for sample in self.metadata.samples:
            if sample.general.bestassemblyfile != 'NA':
                if type(sample[self.analysistype].allelenames) == list:
                    #
                    if sample[self.analysistype].profile != 'NA':
                        # Initialise dictionaries
                        sample[self.analysistype].profilematches = dict()
                        sample[self.analysistype].sequencetypematches = dict()
                        # Create the profiledata variable to avoid writing self.profiledata[self.analysistype]
                        profiledata = sample[self.analysistype].profiledata
                        # For each gene
                        for gene in sorted(sample[self.analysistype].allelenames):
                            try:
                                allelenumber = sample[self.analysistype].allelematches[gene].split('-')[1]
                                # Find the profile with the most alleles in common with the query genome
                                for sequencetype in profiledata:
                                    # refallele is the allele number of the sequence type
                                    refallele = profiledata[sequencetype][gene]
                                    if allelenumber == refallele:
                                        # Add matching alleles
                                        try:
                                            sample[self.analysistype].profilematches[sequencetype] += 1
                                            sample[self.analysistype].sequencetypematches[sequencetype].append(
                                                refallele)
                                        except KeyError:
                                            sample[self.analysistype].profilematches[sequencetype] = 1
                                            sample[self.analysistype].sequencetypematches[sequencetype] = list()
                                            sample[self.analysistype].sequencetypematches[sequencetype].append(
                                                refallele)
                            except KeyError:
                                pass
1242160
def parse_insertion(insertion, gff):
    
    offset = insertion['offset']
    for ins in parse_fasta(insertion['insertion sequence'].split('|')):
        strand = insertion['strand']
        ID = ins[0].split('>')[1].split()[0]
        Start, End = [int(i) for i in ins[0].split('gene-pos=', 1)[1].split()[0].split('-')]
        Start, End = abs(Start + offset), abs(End + offset)
        if strand == '-':
            Start, End = End, Start
        gff['#seqname'].append(insertion['ID'])
        gff['source'].append(insertion['source'])
        gff['feature'].append('IVS')
        gff['start'].append(Start)
        gff['end'].append(End)
        gff['score'].append('.')
        gff['strand'].append(strand) # same as rRNA
        gff['frame'].append('.')
        gff['attribute'].append('ID=%s' % (ID))
    return gff
719558
def get_patient_mhc_haplotype(job, patient_dict):
    
    haplotype_archive = job.fileStore.readGlobalFile(patient_dict['hla_haplotype_files'])
    haplotype_archive = untargz(haplotype_archive, os.getcwd())
    output_dict = {}
    for filename in 'mhci_alleles.list', 'mhcii_alleles.list':
        output_dict[filename] = job.fileStore.writeGlobalFile(os.path.join(haplotype_archive,
                                                                           filename))
    return output_dict
707707
def map_to_precursors_on_fly(seqs, names, loci, args):
    """"""
    precursor = precursor_sequence(loci, args.ref).upper()
    dat = dict()
    for s, n in itertools.izip(seqs, names):
        res = pyMatch.Match(precursor, str(s), 1, 3)
        if res > -1:
            dat[n] = [res, res + len(s)]
    logger.debug("mapped in %s: %s out of %s" % (loci, len(dat), len(seqs)))
    return dat
175037
def make_sequence(seq, name="S"):
    
    return ["{}_{}_{}".format(name, i, x) for i, x in enumerate(seq)]
176247
def load_sequences_to_reference(self, sc=None, force_rerun=False):
        """"""
        log.info('Loading sequences to reference GEM-PRO...')
        from random import shuffle
        g_ids = [g.id for g in self.reference_gempro.functional_genes]
        shuffle(g_ids)

        def _load_sequences_to_reference_gene_sc(g_id, outdir=self.sequences_by_gene_dir,
                                                 g_to_pickle=self.gene_protein_pickles,
                                                 strain_infodict=self.strain_infodict,
                                                 orth_matrix=self.df_orthology_matrix, force_rerun=force_rerun):
            """Load orthologous strain sequences to reference Protein object, save as new pickle"""
            import ssbio.utils
            import ssbio.io
            from Bio import SeqIO
            import os.path as op

            protein_seqs_pickle_path = op.join(outdir, '{}_protein_withseqs.pckl'.format(g_id))

            if ssbio.utils.force_rerun(flag=force_rerun, outfile=protein_seqs_pickle_path):
                protein_pickle_path = g_to_pickle[g_id]
                protein_pickle = ssbio.io.load_pickle(protein_pickle_path)

                for strain, info in strain_infodict.items():
                    strain_sequences = SeqIO.index(info['genome_path'], 'fasta')
                    strain_gene_functional = info['functional_genes'][g_id]
                    if strain_gene_functional:
                        # Pull the gene ID of the strain from the orthology matrix
                        strain_gene_key = orth_matrix.at[g_id, strain]
                        new_id = '{}_{}'.format(g_id, strain)
                        if protein_pickle.sequences.has_id(new_id):
                            continue
                        protein_pickle.load_manual_sequence(seq=strain_sequences[strain_gene_key],
                                                            ident=new_id,
                                                            set_as_representative=False)
                protein_pickle.save_pickle(outfile=protein_seqs_pickle_path)

            return g_id, protein_seqs_pickle_path

        if sc:
            genes_rdd = sc.parallelize(g_ids)
            result = genes_rdd.map(_load_sequences_to_reference_gene_sc).collect()
        else:
            result = []
            for g in tqdm(g_ids):
                result.append(self._load_sequences_to_reference_gene(g, force_rerun))

        log.info('Storing paths to new Protein objects in self.gene_protein_pickles...')
        updated = []
        for g_id, protein_pickle in result:
            self.gene_protein_pickles[g_id] = protein_pickle
            updated.append(g_id)
        not_updated = set(list(self.gene_protein_pickles.keys())).difference(updated)
        log.info('No change to {} genes, removing from gene_protein_pickles'.format(len(not_updated)))
        log.debug(not_updated)
        for rem in not_updated:
            del self.gene_protein_pickles[rem]
172990
def GeneReader( fh, format='gff' ):
    

    known_formats = ( 'gff', 'gtf', 'bed')
    if format not in known_formats: 
        print('%s format not in %s' % (format, ",".join( known_formats )), file=sys.stderr)
        raise Exception('?')

    if format == 'bed':
        for line in fh:    
            f = line.strip().split()
            chrom = f[0]
            chrom_start = int(f[1])
            name = f[4]
            strand = f[5]
            cdsStart = int(f[6])
            cdsEnd = int(f[7])
            blockCount = int(f[9])
            blockSizes = [ int(i) for i in f[10].strip(',').split(',') ]
            blockStarts = [ chrom_start + int(i) for i in f[11].strip(',').split(',') ]

            # grab cdsStart - cdsEnd
            gene_exons = []
            for base,offset in zip( blockStarts, blockSizes ):
                exon_start = base
                exon_end = base+offset
                gene_exons.append( (exon_start, exon_end) )
            
    genelist = {}
    grouplist = []
    if format == 'gff' or format == 'gtf':
        for line in fh:
            if line.startswith('#'): continue
            fields = line.strip().split('\t')
            if len( fields ) < 9: continue

            # fields

            chrom = fields[0]
            ex_st = int( fields[3] ) - 1 # make zero-centered
            ex_end = int( fields[4] ) #+ 1 # make exclusive
            strand = fields[6]

            if format == 'gtf':
                group = fields[8].split(';')[0]
            else:
                group = fields[8]

            if group not in grouplist: grouplist.append( group )
            if group not in genelist:
                genelist[group] = (chrom, strand, [])
            exons_i = 2
            genelist[group][exons_i].append( ( ex_st, ex_end ) )

        sp = lambda a,b: cmp( a[0], b[0] )

        #for gene in genelist.values():
        for gene in grouplist:
            chrom, strand, gene_exons = genelist[ gene ]
            gene_exons = bitset_union( gene_exons )
            yield chrom, strand, gene_exons, gene
406916
def create_fa(self):
        

        if self._seqs is None:
            os.symlink(self._fa0_fn, self._fa_fn)
        else:
            in_seqs = pyfaidx.Fasta(self._fa0_fn)

            with open(self._fa_fn, "w+") as g:
                for seq_desc in self._seqs:
                    x = in_seqs[seq_desc]
                    name, seq = x.name, str(x)
                    g.write(">" + name + "\n")
                    n = 80
                    seq_split = "\n".join([seq[i:i + n] for i in range(0, len(seq), n)])
                    g.write(seq_split + "\n")