Acquire the location of a gene from SGD http://www.yeastgenome.org
    :param gene_name: Name of the gene.
    :type gene_name: string
    :returns location: [int: chromosome, int:biostart, int:bioend, int:strand]
    :rtype location: list
550438
def get_yeast_gene_location(gene_name):
    '''

    '''
    from intermine.webservice import Service
    service = Service('http://yeastmine.yeastgenome.org/yeastmine/service')

    # Get a new query on the class (table) you will be querying:
    query = service.new_query('Gene')

    # The view specifies the output columns
    query.add_view('primaryIdentifier', 'secondaryIdentifier', 'symbol',
                   'name', 'organism.shortName',
                   'chromosome.primaryIdentifier',
                   'chromosomeLocation.start', 'chromosomeLocation.end',
                   'chromosomeLocation.strand')

    # Uncomment and edit the line below (the default) to select a custom sort
    # order:
    # query.add_sort_order('Gene.primaryIdentifier', 'ASC')

    # You can edit the constraint values below
    query.add_constraint('organism.shortName', '=', 'S. cerevisiae',
                         code='B')
    query.add_constraint('Gene', 'LOOKUP', gene_name, code='A')

    # Uncomment and edit the code below to specify your own custom logic:
    # query.set_logic('A and B')
    chromosomes = {'chrI': 1,
                   'chrII': 2,
                   'chrIII': 3,
                   'chrIV': 4,
                   'chrV': 5,
                   'chrVI': 6,
                   'chrVII': 7,
                   'chrVIII': 8,
                   'chrIX': 9,
                   'chrX': 10,
                   'chrXI': 11,
                   'chrXII': 12,
                   'chrXIII': 13,
                   'chrXIV': 14,
                   'chrXV': 15,
                   'chrXVI': 16}
    first_result = query.rows().next()

    return [chromosomes[first_result['chromosome.primaryIdentifier']],
            first_result['chromosomeLocation.start'],
            first_result['chromosomeLocation.end'],
            int(first_result['chromosomeLocation.strand'])]
760098
def gene_names(self):
        
        return self.ensembl.gene_names_at_locus(
            self.contig, self.start, self.end)
976365
def get_transcript_ids_for_ensembl_gene_ids(self, gene_ids, hgnc_symbols):
        
        
        chroms = {"1", "2", "3", "4", "5", "6", "7", "8", "9", "10", "11", \
             "12", "13", "14", "15", "16", "17", "18", "19", "20", "21", "22", \
              "X", "Y"}
        
        headers = {"content-type": "application/json"}
        
        transcript_ids = []
        for gene_id in gene_ids:
            self.attempt = 0
            ext = "/overlap/id/{}?feature=transcript".format(gene_id)
            r = self.ensembl_request(ext, headers)
            
            for item in json.loads(r):
                # ignore non-coding transcripts
                if item["biotype"] not in ["protein_coding", "polymorphic_pseudogene"]:
                    continue
                
                # ignore transcripts not on the standard chromosomes
                # (non-default chroms fail to map the known de novo variants
                # to the gene location
                if item["Parent"] != gene_id or item["seq_region_name"] not in \
                        chroms or \
                        all([symbol not in item["external_name"] for symbol in hgnc_symbols]):
                    continue
                transcript_ids.append(item["id"])
        
        return transcript_ids
809486
def _get_genes(self, variant):
        
        transcripts = variant['transcripts']
        ensembl_ids = [transcript['ensembl_id'] for transcript in
                       transcripts if transcript['ensembl_id']]
        hgnc_symbols = [transcript['hgnc_symbol'] for transcript in
                        transcripts if transcript['hgnc_symbol']]
        genes = get_gene_info(ensembl_ids, hgnc_symbols)
        return genes
356759
def get_genomic_seq_for_region(self, chrom, start_pos, end_pos):
        
        
        headers = {"content-type": "text/plain"}
        
        self.attempt = 0
        ext = "/sequence/region/human/{}:{}..{}:1".format(chrom, start_pos, end_pos)
        
        return self.ensembl_request(ext, headers)
113869
def _load_knownGene(filename):
        
        try:
            df = pd.read_table(
                filename,
                names=[
                    "name",
                    "chrom",
                    "strand",
                    "txStart",
                    "txEnd",
                    "cdsStart",
                    "cdsEnd",
                    "exonCount",
                    "exonStarts",
                    "exonEnds",
                    "proteinID",
                    "alignID",
                ],
                index_col=0,
            )
            df["chrom"] = df["chrom"].str[3:]
            return df
        except Exception as err:
            print(err)
            return None
1326099
def find_gene(self, name: str):
        """"""
        result = [ExpGene.from_series(s)
                  for i, s in self.loc[self['name'] == name].iterrows()]
        return result
31071
def from_seqfeature(s, **kwargs):
    
    source = s.qualifiers.get('source', '.')[0]
    score = s.qualifiers.get('score', '.')[0]
    seqid = s.qualifiers.get('seqid', '.')[0]
    frame = s.qualifiers.get('frame', '.')[0]
    strand = _feature_strand[s.strand]

    # BioPython parses 1-based GenBank positions into 0-based for use within
    # Python.  We need to convert back to 1-based GFF format here.
    start = s.location.start.position + 1
    stop = s.location.end.position
    featuretype = s.type
    id = s.id
    attributes = dict(s.qualifiers)
    attributes.pop('source', '.')
    attributes.pop('score', '.')
    attributes.pop('seqid', '.')
    attributes.pop('frame', '.')
    return Feature(seqid, source, featuretype, start, stop, score, strand,
                   frame, attributes, id=id, **kwargs)
750108
def gene_features(self):
        
        nm, strand = self.gene_name, self.strand
        feats = [(self.chrom, self.start, self.end, nm, strand, 'gene')]
        for feat in ('introns', 'exons', 'utr5', 'utr3', 'cdss'):
            fname = feat[:-1] if feat[-1] == 's' else feat
            res = getattr(self, feat)
            if res is None or all(r is None for r in res): continue
            if not isinstance(res, list): res = [res]
            feats.extend((self.chrom, s, e, nm, strand, fname) for s, e in res)

        tss = self.tss(down=1)
        if tss is not None:
            feats.append((self.chrom, tss[0], tss[1], nm, strand, 'tss'))
            prom = self.promoter()
            feats.append((self.chrom, prom[0], prom[1], nm, strand, 'promoter'))

        return sorted(feats, key=itemgetter(1))
1047823
def get_linc_rna_genes(
        path_or_buffer,
        remove_duplicates=True,
        **kwargs):
    

    """
    valid_biotypes = set(['lincRNA'])
    
    df = get_genes(path_or_buffer, valid_biotypes,
                   remove_duplicates=remove_duplicates, **kwargs)
    return df
254543
def get_coords(data):
    
    for category, vtypes in [("LOH", {"LOSS", "HETEROZYGOSITY"}),
                             ("amplification", {"AMPLIFICATION"})]:
        out = tz.get_in([category, dd.get_genome_build(data)], _COORDS, {})
        priority_file = dd.get_svprioritize(data)
        if priority_file:
            if os.path.basename(priority_file).find("civic") >= 0:
                for chrom, start, end, gene in _civic_regions(priority_file, vtypes, dd.get_disease(data)):
                    out[gene] = (chrom, start, end)
            elif os.path.basename(priority_file).find(".bed") >= 0:
                for line in utils.open_gzipsafe(priority_file):
                    parts = line.strip().split("\t")
                    if len(parts) >= 4:
                        chrom, start, end, gene = parts[:4]
                        out[gene] = (chrom, int(start), int(end))
        yield category, out
1257589
def get_gene2section2gos(gene2gos, sec2gos):
        """"""
        gene2section2gos = {}
        for geneid, gos_gene in gene2gos.items():
            section2gos = {}
            for section_name, gos_sec in sec2gos.items():
                gos_secgene = gos_gene.intersection(gos_sec)
                if gos_secgene:
                    section2gos[section_name] = gos_secgene
            gene2section2gos[geneid] = section2gos
        return gene2section2gos
980353
def set_gene_name(self,name):
    
    self._options = self._options._replace(gene_name = name)
321299
def init_genome_coordinates(self) :
        """"""
        self.seqpos2genome = {}

        # record genome positions for each sequence position
        seq_pos = 0
        for estart, eend in self.exons:
            for genome_pos in range(estart, eend):
                if self.strand == '+':
                    self.seqpos2genome[seq_pos] = genome_pos
                elif self.strand == '-':
                    tmp = self.cds_len - seq_pos - 1
                    self.seqpos2genome[tmp] = genome_pos
                seq_pos += 1

        # recode 5' splice site locations
        for i in range(0, self.five_ss_len):
            seq_pos = self.cds_len + i
            ss_ix = i // 2  # the ss_ix'th 5'ss starting from upstream tx
            pos_in_ss = i % 2  # whether first/second nuc in splice site

            # determine genome coordinates for 5' splice site
            if self.strand == '+':
                self.seqpos2genome[seq_pos] = self.exons[ss_ix][1] + pos_in_ss
            else:
                exon_pos = -1 - ss_ix
                self.seqpos2genome[seq_pos] = self.exons[exon_pos][0] - pos_in_ss - 1

        # recode 3' splice site locations
        for i in range(0, self.three_ss_len):
            seq_pos = self.cds_len + self.five_ss_len + i
            ss_ix = i // 2  # the ss_ix'th 3'ss starting from upstream tx
            pos_in_ss = i % 2  # whether first/second nuc in splice site

            # determine genome coordinates for 3' splice site
            if self.strand == '+':
                self.seqpos2genome[seq_pos] = self.exons[ss_ix+1][0] - 2 + pos_in_ss
            else:
                exon_pos = -1 - ss_ix
                self.seqpos2genome[seq_pos] = self.exons[exon_pos-1][1] + 1 - pos_in_ss
347100
def _biotype_lookup_fn(gtf):
    
    db = get_gtf_db(gtf)
    sources = set([feature.source for feature in db.all_features()])
    gene_biotypes = set([feature.attributes.get("gene_biotype", [None])[0]
                         for feature in db.all_features()])
    biotypes = set([feature.attributes.get("biotype", [None])[0]
                    for feature in db.all_features()])
    if "protein_coding" in sources:
        return lambda feature: feature.source
    elif "protein_coding" in biotypes:
        return lambda feature: feature.attributes.get("biotype", [None])[0]
    elif "protein_coding" in gene_biotypes:
        return lambda feature: feature.attributes.get("gene_biotype", [None])[0]
    else:
        return None
1180194
def get_gene2aart(gene2section2gos, sec2chr):
        """"""
        geneid2str = {}
        for geneid, section2gos_gene in gene2section2gos.items():
            letters = [abc if s in section2gos_gene else "." for s, abc in sec2chr.items()]
            geneid2str[geneid] = "".join(letters)
        return geneid2str
470748
def referencegenomefinder(self):
        
        # Initialise dictionaries
        referencematch = defaultdict(make_dict)
        referencehits = defaultdict(make_dict)
        # Set the name of the reference profile file
        referencegenomeprofile = '{}rMLST_referenceprofile.json'.format(self.referenceprofilepath)
        # Open the reference profile and load the profile into memory
        with open(referencegenomeprofile) as referencefile:
            referencetypes = json.load(referencefile)
        # Iterate through the samples
        for sample in self.metadata:
            if sample[self.analysistype].reportdir != 'NA':
                # Iterate through the reference genomes in the profile
                for genome in referencetypes:
                    # Initialise the number of identical alleles between the assembly of interest and the
                    # reference genome to 0
                    referencehits[sample.name][genome] = 0
                    # Iterate through all the genes in the analysis
                    for gene in self.bestdict[sample.name]:
                        # If the alleles match between the assembly of interest and the reference genome, increment
                        # the number of matches by 1
                        if list(self.bestdict[sample.name][gene].keys())[0] == referencetypes[genome][gene]:
                            referencematch[sample.name][genome][gene] = 1
                            referencehits[sample.name][genome] += 1
                        else:
                            referencematch[sample.name][genome][gene] = 0
        #
        for sample in self.metadata:
            if sample[self.analysistype].reportdir != 'NA':
                # Get the best number of matches
                # From: https://stackoverflow.com/questions/613183/sort-a-python-dictionary-by-value
                try:
                    matches = sorted(referencehits[sample.name].items(),
                                     key=operator.itemgetter(1), reverse=True)
                    most_matches = matches[0][1]
                    i = 0
                    match_list = list()
                    while matches[i][1] == most_matches:
                        match_list.append(matches[i])
                        i += 1
                    sorted_list = sorted(match_list)
                    sortedmatches = sorted_list[0]
                except IndexError:
                    sortedmatches = (0, 0)
                # If there are fewer matches than the total number of genes in the typing scheme
                if 0 < int(sortedmatches[1]) < len(sample[self.analysistype].allelenames):
                    mismatches = []
                    # Iterate through the gene in the analysis
                    for gene, allele in referencetypes[sortedmatches[0]].items():
                        # Populate :self.referencegenome with the genome name, best reference match, number of matches,
                        # gene, query allele(s), and percent identity
                        percentidentity = '{:.2f}'.format(list(self.bestdict[sample.name][gene].values())[0])
                        self.referencegenome[sample.name][sortedmatches[0]][sortedmatches[1]][gene][list(self.bestdict[
                            sample.name][gene].keys())[0]] = percentidentity
                        if list(self.bestdict[sample.name][gene].keys())[0] != allele:
                            sample[self.analysistype].referencegenome = sortedmatches[0]
                            sample.general.referencegenus = sortedmatches[0].split('_')[0]
                            sample[self.analysistype].referencegenomepath = '{}{}.fa' \
                                .format(self.referenceprofilepath, sortedmatches[0])
                            sample[self.analysistype].matchestoreferencegenome = sortedmatches[1]
                            mismatches.append(({gene: ('{} ({})'.format(list(self.bestdict[sample.name][gene]
                                                                        .keys())[0], allele))}))
                        # sample[self.analysistype].mismatchestoreferencegenome = sorted(mismatches)
                        sample[self.analysistype].mismatchestoreferencegenome = mismatches
                elif sortedmatches == 0:
                    for gene in sample[self.analysistype].allelenames:
                        # Populate the profile of results with 'negative' values for sequence type and sorted matches
                        self.referencegenome[sample.name][sortedmatches[0]][0][gene]['NA'] = 0
                        sample[self.analysistype].referencegenome = 'NA'
                        sample.general.referencegenus = 'NA'
                        sample[self.analysistype].referencegenomepath = 'NA'
                        sample[self.analysistype].matchestoreferencegenome = 0
                        sample[self.analysistype].mismatchestoreferencegenome = [0]
                # Otherwise, the query profile matches the reference profile
                else:
                    for gene in referencetypes[sortedmatches[0]]:
                        # Populate self.referencegenome as above
                        self.referencegenome[sample.name][sortedmatches[0]][sortedmatches[1]][gene][list(self.bestdict[
                            sample.name][gene].keys())[0]] = '{:.2f}'.format(list(self.bestdict[
                                                                                sample.name][gene].values())[0])
                        sample[self.analysistype].referencegenome = sortedmatches[0]
                        sample[self.analysistype].referencegenomepath = '{}{}.fa' \
                            .format(self.referenceprofilepath, sortedmatches[0])
                        sample.general.referencegenus = sortedmatches[0].split('_')[0]
                        sample[self.analysistype].matchestoreferencegenome = sortedmatches[1]
                        sample[self.analysistype].mismatchestoreferencegenome = [0]
        # Print the results to file
        make_path(self.reportpath)
        with open('{}referencegenomes.csv'.format(self.reportpath), 'w') as referencegenomereport:
            row = 'Strain,referencegenome,numberofmatches\n'
            for sample in self.metadata:
                if sample[self.analysistype].reportdir != 'NA':
                    row += '{},{},{}\n'.format(sample.name, sample[self.analysistype].referencegenome,
                                               sample[self.analysistype].matchestoreferencegenome)
            referencegenomereport.write(row)
            dotter()