list: List of genes with no mapping to UniProt.
426084
def missing_uniprot_mapping(self):
        """"""
        uniprot_missing = []
        for g in self.genes:
            ups = g.protein.filter_sequences(UniProtProp)

            no_sequence_file_available = True
            for u in ups:
                if u.sequence_file or u.metadata_file:
                    no_sequence_file_available = False
                    break
            if no_sequence_file_available:
                uniprot_missing.append(g.id)

        return list(set(uniprot_missing))
876452
def find_features(seqs, locus_tag="all", utr_len=200):
    """"""
    found_features = []

    for seq_i in seqs:
        for feature in seq_i.features:
            if feature.type == "CDS" and (locus_tag == "all" or \
                    ('locus_tag' in feature.qualifiers and \
                    feature.qualifiers['locus_tag'][0] == locus_tag)):
                start = max(0, feature.location.nofuzzy_start - utr_len)
                stop  = max(0, feature.location.nofuzzy_end + utr_len)
                feature_seq = seq_i.seq[start:stop]
                f_match = FeatureMatch(feature, feature_seq, feature.strand,
                                       utr_len)
                found_features.append(f_match)

    return found_features
1271376
def convert_uniprot_to_entrez(self, uniprot):
        """"""
        # Submit request to NCBI eutils/Gene Database
        server = "http://eutils.ncbi.nlm.nih.gov/entrez/eutils/esearch.fcgi?" + self.options + "&db=gene&term={0}".format(
            uniprot)
        r = requests.get(server, headers={"Content-Type": "text/xml"})
        if not r.ok:
            r.raise_for_status()
            sys.exit()
        # Process Request
        response = r.text
        info = xmltodict.parse(response)
        geneId = info['eSearchResult']['IdList']['Id']
        # check to see if more than one result is returned
        # if you have more than more result then check which Entrez Id returns the same uniprot Id entered.
        if len(geneId) > 1:
            for x in geneId:
                c = self.convert_entrez_to_uniprot(x)
                c = c.lower()
                u = uniprot.lower()
                if c == u:
                    return x
        else:
            return geneId
1119459
def list_available_genomes(provider=None):
    
    if provider:
        providers = [ProviderBase.create(provider)]
    else:
        # if provider is not specified search all providers
        providers = [ProviderBase.create(p) for 
                        p in ProviderBase.list_providers()]

    for p in providers:
        for row in p.list_available_genomes():
            yield [p.name] + list(row)
144949
def df_representative_sequences(self):
        """"""

        seq_mapping_pre_df = []
        df_cols = ['gene', 'uniprot', 'kegg', 'num_pdbs', 'pdbs', 'seq_len', 'sequence_file', 'metadata_file']

        for g in self.genes_with_a_representative_sequence:
            gene_dict = g.protein.representative_sequence.get_dict(df_format=True, only_attributes=df_cols)
            gene_dict['gene'] = g.id
            seq_mapping_pre_df.append(gene_dict)

        df = pd.DataFrame.from_records(seq_mapping_pre_df, columns=df_cols).set_index('gene')
        if df.empty:
            log.warning('Empty dataframe')
            return df
        else:
            return ssbio.utils.clean_df(df)
683988
def gffatributes2ids(s):
    
    Name,gene_id,transcript_id,protein_id,exon_id=np.nan,np.nan,np.nan,np.nan,np.nan
    if '=' in s:
        d=dict([i.split('=') for i in s.split(';')])
        if 'Parent' in d:
            d[d['Parent'].split(':')[0]+'_id']=d['Parent'].split(':')[1]
        Name,gene_id,transcript_id,protein_id,exon_id=np.nan,np.nan,np.nan,np.nan,np.nan
        if 'Name' in d:    
            Name=d['Name']
        if 'gene_id' in d:    
            gene_id=d['gene_id']
        if 'transcript_id' in d:    
            transcript_id=d['transcript_id']
        if 'protein_id' in d:    
            protein_id=d['protein_id']
        if 'exon_id' in d:    
            exon_id=d['exon_id']
    return Name,gene_id,transcript_id,protein_id,exon_id
1053895
def dalignbedannot2daggbyquery(cfg):
    
    datatmpd=cfg['datatmpd']

    dalignbedannot=del_Unnamed(pd.read_csv(cfg['dalignbedannotp'],sep='\t',low_memory=False))
    
    daggbyqueryp=f'{datatmpd}/10_daggbyquery.tsv'      
    logging.info(basename(daggbyqueryp))
    if not exists(daggbyqueryp) or cfg['force']:
        dalignbedannot=dfliststr2dflist(dalignbedannot,
                                ['type', 'gene name', 'gene id', 'transcript id', 'protein id', 'exon id'],
                                colfmt='tuple')    
        dalignbedannot['alternate alignments count']=1
        import itertools
        from rohan.dandage.io_sets import unique
        def unique_dropna(l): return unique(l,drop='nan')
        def merge_unique_dropna(l): return unique(list(itertools.chain(*l)),drop='nan')
        cols2aggf={'id':unique_dropna,
         'type':merge_unique_dropna,
         'gene name':merge_unique_dropna,
         'gene id':merge_unique_dropna,
         'transcript id':merge_unique_dropna,
         'protein id':merge_unique_dropna,
         'exon id':merge_unique_dropna,
         'region':unique_dropna,
         'alternate alignments count':sum,
        }
        daggbyquery=dalignbedannot.groupby('query id').agg(cols2aggf)
        daggbyquery.to_csv(daggbyqueryp,sep='\t')
        daggbyquery.to_csv(cfg['dalignannotedp'],sep='\t')
    return cfg
1003262
def fetch_protein_gene_map(self, taxon_id):
        
        protein_dict = dict()
        # col = self.columns['ensembl_biomart']
        col = ['ensembl_peptide_id', 'ensembl_gene_id']
        raw_query = self._build_biomart_gene_query(taxon_id, col)

        params = urllib.parse.urlencode({'query': raw_query})

        conn = http.client.HTTPConnection(ENS_URL)
        conn.request("GET", '/biomart/martservice?' + params)
        response = conn.getresponse()

        for line in response:
            line = line.decode('utf-8').rstrip()
            row = line.split('\t')

            if len(row) != len(col):
                LOG.warning("Data error for p2g query on %s", taxon_id)
                LOG.warning("Expected columns for \n%s", col)
                LOG.warning("Got data \n%s", row)
                continue

            protein_dict[
                row[col.index('ensembl_peptide_id')]] = [
                    row[col.index('ensembl_gene_id')]
                ]

        conn.close()
        LOG.info(
            "length protien list for taxon: %s is %i", taxon_id, len(protein_dict))
        return protein_dict
1057932
def print_table(seqs, id2name, name):
    
    itable = open('%s.itable' % (name.rsplit('.', 1)[0]), 'w')
    print('\t'.join(['#sequence', 'gene', 'model', 'insertion', 'gene position', 'model position', 'length', 'orf?', 'intron?', 'orf?intron?', 'insertion', 'orf', 'intron']), file=itable)
    for seq, info in list(seqs.items()):
        gene, model, insertions = info
        name = id2name[seq]
        for i, ins in enumerate(insertions, 1):
            gene_pos, model_pos, length, iseq, \
                    orfs, introns, orfs_b, introns_b, orf_annotations = ins
            # append annotation to orf header
            for orf in orfs:
                parts = orf[0].split()
                annotation = orf_annotations[parts[0].split('>')[1]]
                orf[0] = '%s %s %s' % (parts[0], annotation, ' '.join(parts[1:]))
            # get orf position
            gene_pos = '-'.join([str(j) for j in gene_pos])
            # check if orf, intron is present
            if orfs_b is True or introns_b is True:
                orfs_introns_b = True
            else:
                orfs_introns_b = False
            out = [name, gene, model, i, gene_pos, model_pos, length, orfs_b, introns_b, orfs_introns_b]
            out.append('|'.join(iseq))
            out.append('|'.join(['|'.join(orf) for orf in orfs]))
            out.append('|'.join(['|'.join(intron) for intron in introns]))
            print('\t'.join([str(i) for i in out]), file=itable)
    itable.close()
461018
def get_progenies(p1, p2, x_linked=False, tolerance=0):
    
    _p1 = expand_alleles(p1, tolerance=tolerance)
    _p2 = expand_alleles(p2, tolerance=tolerance)
    possible_progenies = set(tuple(sorted(x)) for x in product(_p1, _p2))
    if x_linked:  # Add all hemizygotes
        possible_progenies |= set((x, x) for x in (set(_p1) | set(_p2)))
    return possible_progenies
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
237241
def blat(self, db=None, sequence=None, seq_type="DNA"):
        
        from . blat_blast import blat, blat_all
        assert sequence in (None, "cds", "mrna")
        seq = self.sequence() if sequence is None else ("".join(self.cds_sequence if sequence == "cds" else self.mrna_sequence))
        if isinstance(db, (tuple, list)):
            return blat_all(seq, self.gene_name, db, seq_type)
        else:
            return blat(seq, self.gene_name, db or self.db, seq_type)
155014
static public void addGmodGFF3GeneFeatures(LinkedHashMap<String, ChromosomeSequence> chromosomeSequenceList, FeatureList listGenes) throws Exception {


		// key off mRNA as being a known feature that may or may not have a parent gene


		FeatureList mRNAFeatures = listGenes.selectByType("mRNA");
		LinkedHashMap<String,FeatureList> featureIDHashMap = FeatureHelper.buildFeatureAtrributeIndex("ID", listGenes);
		LinkedHashMap<String,FeatureList> featureParentHashMap = FeatureHelper.buildFeatureAtrributeIndex("Parent", listGenes);

		for (FeatureI f : mRNAFeatures) {
			String geneID;
			String geneNote = null;
			String geneSource = null;
			String sequenceName = null;
			ChromosomeSequence seq = null;
			GeneSequence geneSequence = null;

			Feature mRNAFeature = (Feature) f;
			String mRNAID = mRNAFeature.getAttribute("ID");
			String mRNAsource = mRNAFeature.source();
			String mRNANote = mRNAFeature.getAttribute("Note");
			String mRNAParent = mRNAFeature.getAttribute("Parent");
			if (mRNAParent != null && mRNAParent.length() > 0) {
			   // FeatureList geneFeatureList = listGenes.selectByAttribute("ID", mRNAParent);
				FeatureList geneFeatureList = featureIDHashMap.get(mRNAParent);
				Feature geneFeature = (Feature) geneFeatureList.get(0);
				geneID = geneFeature.getAttribute("ID");
				geneNote = geneFeature.getAttribute("Note");
				geneSource = geneFeature.source();
				sequenceName = geneFeature.seqname();

				//
			} else {
				//deal with cases where no parent gene is given
				geneID = mRNAID;
				geneSource = mRNAsource;
				sequenceName = mRNAFeature.seqname();
			}

			seq = chromosomeSequenceList.get(sequenceName);

			AccessionID geneAccessionID = new AccessionID(geneID);

		  //  FeatureList mRNAChildren = listGenes.selectByAttribute("Parent", mRNAID);
			FeatureList mRNAChildren = featureParentHashMap.get(mRNAID);
			FeatureList cdsFeatures = mRNAChildren.selectByType("CDS");
			FeatureI feature = cdsFeatures.get(0);
			Strand strand = Strand.POSITIVE;

			if (feature.location().isNegative()) {
				strand = Strand.NEGATIVE;
			}
			cdsFeatures = cdsFeatures.sortByStart();







			//String seqName = feature.seqname();
			FeatureI startCodon = null;
			FeatureI stopCodon = null;
			Integer startCodonBegin = null;
			Integer stopCodonEnd = null;
			String startCodonName = "";
			String stopCodonName = "";
			FeatureList startCodonList = mRNAChildren.selectByType("five_prime_UTR");
			if (startCodonList != null && startCodonList.size() > 0) {
				startCodon = startCodonList.get(0);
				if (strand == Strand.NEGATIVE) {
					startCodonBegin = startCodon.location().bioEnd();
				} else {
					startCodonBegin = startCodon.location().bioStart();
				}
				startCodonName = startCodon.getAttribute("ID");
			}

			FeatureList stopCodonList = mRNAChildren.selectByType("three_prime_UTR");

			if (stopCodonList != null && stopCodonList.size() > 0) {
				stopCodon = stopCodonList.get(0);
				if (strand == Strand.NEGATIVE) {
					stopCodonEnd = stopCodon.location().bioStart();
				} else {
					stopCodonEnd = stopCodon.location().bioEnd();
				}
				stopCodonName = stopCodon.getAttribute("ID");

			}




			if (startCodonBegin == null) {
				if (strand == Strand.NEGATIVE) {
					FeatureI firstFeature = cdsFeatures.get(0);

					startCodonBegin = firstFeature.location().bioEnd();
				} else {
					FeatureI firstFeature = cdsFeatures.get(0);

					startCodonBegin = firstFeature.location().bioStart();
				}
			}

			if (stopCodonEnd == null) {
				if (strand == Strand.NEGATIVE) {
					FeatureI lastFeature = cdsFeatures.get(cdsFeatures.size() - 1);
					stopCodonEnd = lastFeature.location().bioStart();
				} else {
					FeatureI lastFeature = cdsFeatures.get(cdsFeatures.size() - 1);
					stopCodonEnd = lastFeature.location().bioEnd();
				}
			}
			//for gtf ordering can be strand based so first is last and last is first
			if (startCodonBegin > stopCodonEnd) {
				int temp = startCodonBegin;
				startCodonBegin = stopCodonEnd;
				stopCodonEnd = temp;
			}



			AccessionID transcriptAccessionID = new AccessionID(mRNAID);
			geneSequence = seq.getGene(geneID);
			if (geneSequence == null) {
				geneSequence = seq.addGene(geneAccessionID, startCodonBegin, stopCodonEnd, strand);
				geneSequence.setSource(geneSource);
				if (geneNote != null && geneNote.length() > 0) {
					geneSequence.addNote(geneNote);
				}
			} else {

				if (startCodonBegin < geneSequence.getBioBegin()) {
					geneSequence.setBioBegin(startCodonBegin);
				}
				if (stopCodonEnd > geneSequence.getBioBegin()) {
					geneSequence.setBioEnd(stopCodonEnd);
				}

			}
			TranscriptSequence transcriptSequence = geneSequence.addTranscript(transcriptAccessionID, startCodonBegin, stopCodonEnd);
			transcriptSequence.setSource(mRNAsource);
			if (mRNANote != null && mRNANote.length() > 0) {
				transcriptSequence.addNote(mRNANote);

			}
			if (startCodon != null) {
				if (startCodonName == null || startCodonName.length() == 0) {
					startCodonName = mRNAID + "-start_codon-" + startCodon.location().bioStart() + "-" + startCodon.location().bioEnd();
				}
				transcriptSequence.addStartCodonSequence(new AccessionID(startCodonName), startCodon.location().bioStart(), startCodon.location().bioEnd());
			}
			if (stopCodon != null) {
				if (stopCodonName == null || stopCodonName.length() == 0) {
					stopCodonName = mRNAID + "-stop_codon-" + stopCodon.location().bioStart() + "-" + stopCodon.location().bioEnd();
				}
				transcriptSequence.addStopCodonSequence(new AccessionID(stopCodonName), stopCodon.location().bioStart(), stopCodon.location().bioEnd());
			}

			for (FeatureI cdsFeature : cdsFeatures) {
				Feature cds = (Feature) cdsFeature;
				String cdsNote = cdsFeature.getAttribute("Note");
				String cdsSource = cds.source();
				String cdsName = cds.getAttribute("ID");
				if (cdsName == null || cdsName.length() == 0) {
					cdsName = mRNAID + "-cds-" + cds.location().bioStart() + "-" + cds.location().bioEnd();
				}
				AccessionID cdsAccessionID = new AccessionID(cdsName);
				ExonSequence exonSequence = geneSequence.addExon(cdsAccessionID, cdsFeature.location().bioStart(), cdsFeature.location().bioEnd());
				exonSequence.setSource(cdsSource);
				if (cdsNote != null && cdsNote.length() > 0) {
					exonSequence.addNote(cdsNote);
				}
				transcriptSequence.addCDS(cdsAccessionID, cdsFeature.location().bioStart(), cdsFeature.location().bioEnd(), cds.frame());
			}
			geneSequence.addIntronsUsingExons();

		}

	}
1295844
def parse_star_genecount_report(self, f):
        
        # Three numeric columns: unstranded, stranded/first-strand, stranded/second-strand
        keys = [ 'N_unmapped', 'N_multimapping', 'N_noFeature', 'N_ambiguous' ]
        unstranded = { 'N_genes': 0 }
        first_strand = { 'N_genes': 0 }
        second_strand = { 'N_genes': 0 }
        num_errors = 0
        num_genes = 0
        for l in f['f']:
            s = l.split("\t")
            try:
                for i in [1,2,3]:
                    s[i] = float(s[i])
                if s[0] in keys:
                    unstranded[s[0]] = s[1]
                    first_strand[s[0]] = s[2]
                    second_strand[s[0]] = s[3]
                else:
                    unstranded['N_genes'] += s[1]
                    first_strand['N_genes'] += s[2]
                    second_strand['N_genes'] += s[3]
                    num_genes += 1
            except IndexError:
                # Tolerate a few errors in case there is something random added at the top of the file
                num_errors += 1
                if num_errors > 10 and num_genes == 0:
                    log.warning("Error parsing {}".format(f['fn']))
                    return None
        if num_genes > 0:
            return { 'unstranded': unstranded, 'first_strand': first_strand, 'second_strand': second_strand }
        else:
            return None
270281
def getGenomeList() :
	""""""
	import rabaDB.filters as rfilt
	f = rfilt.RabaQuery(Genome_Raba)
	names = []
	for g in f.iterRun() :
		names.append(g.name)
	return names
911466
@Override
	public void parseHeader(String header, S sequence) {
		//uniptrot
		// tr|Q0TET7|Q0TET7_ECOL5 Putative uncharacterized protein OS=Escherichia coli O6:K15:H31 (strain 536 / UPEC) GN=ECP_2553 PE=4 SV=1
		sequence.setOriginalHeader(header);
		String[] data = getHeaderValues(header);

		if (data.length == 1) {
			sequence.setAccession(new AccessionID(data[0]));
		} else  if (data[0].equalsIgnoreCase("sp") || data[0].equalsIgnoreCase("tr")) {
			if (data[0].equalsIgnoreCase("sp")) {
				sequence.setAnnotationType(AnnotationType.CURATED);
			} else {
				sequence.setAnnotationType(AnnotationType.PREDICTED);
			}

			sequence.setAccession(new AccessionID(data[1], DataSource.UNIPROT));
			if (data.length > 2) {
				sequence.setDescription(data[2]);
			}

		} else if (data[0].equalsIgnoreCase("gi")) {
			DataSource giSource = DataSource.UNKNOWN;
			if (data.length >= 3) {
				if (data[2].equalsIgnoreCase("gb")) {
					giSource = DataSource.GENBANK;
				} else if (data[2].equalsIgnoreCase("emb")) {
					giSource = DataSource.ENA;
				} else if (data[2].equalsIgnoreCase("dbj")) {
					giSource = DataSource.DDBJ;
				}
				sequence.setAccession(new AccessionID(data[3], giSource));
			} else {
				sequence.setAccession(new AccessionID(header, giSource));
			}
		} else if (data[0].equalsIgnoreCase("pir")) {
			sequence.setAccession(new AccessionID(data[2], DataSource.NBRF));
		} else if (data[0].equalsIgnoreCase("prf")) {
			sequence.setAccession(new AccessionID(data[2], DataSource.PRF));
		} else if (data[0].equalsIgnoreCase("pdb")) {
			sequence.setAccession(new AccessionID(data[1] + ":" + data[2], DataSource.PDB1));
		} else if (data[0].startsWith("PDB")) {
			String[] pdbe = data[0].split(" ");
			String[] pdbaccession = pdbe[0].split(":");
			sequence.setAccession(new AccessionID(pdbaccession[1], DataSource.PDBe));
		} else if (data[0].indexOf(":") != -1 && data.length > 1 && data[1].equals("PDBID")) {
			sequence.setAccession(new AccessionID(data[0], DataSource.PDB2));
		} else if (data[0].equalsIgnoreCase("pat")) {
			sequence.setAccession(new AccessionID(data[2], DataSource.PATENTS));
		} else if (data[0].equalsIgnoreCase("bbs")) {
			sequence.setAccession(new AccessionID(data[1], DataSource.GENINFO));
		} else if (data[0].equalsIgnoreCase("gnl")) {
			sequence.setAccession(new AccessionID(data[2], DataSource.GENERAL));
		} else if (data[0].equalsIgnoreCase("ref")) {
			sequence.setAccession(new AccessionID(data[1], DataSource.NCBI));
		} else if (data[0].equalsIgnoreCase("lcl")) {
			sequence.setAccession(new AccessionID(data[1], DataSource.LOCAL));
		} else {
			sequence.setAccession(new AccessionID(data[0])); // avoid the common problem of picking up all the comments original header in getOriginalHeader
		}


	}
57053
def _align_orthologous_gene_pairwise(self, g_id, gapopen=10, gapextend=0.5, engine='needle', parse=True, force_rerun=False):
        """"""
        protein_seqs_aln_pickle_path = op.join(self.sequences_by_gene_dir, '{}_protein_withseqs_dis_aln.pckl'.format(g_id))

        if ssbio.utils.force_rerun(flag=force_rerun, outfile=protein_seqs_aln_pickle_path):
            protein_seqs_pickle_path = self.gene_protein_pickles[g_id]
            protein_pickle = ssbio.io.load_pickle(protein_seqs_pickle_path)

            if not protein_pickle.representative_sequence:
                log.error('{}: no representative sequence to align to'.format(g_id))
                return

            if len(protein_pickle.sequences) < 1:
                log.error('{}: no other sequences to align to'.format(g_id))
                return

            alignment_dir = op.join(self.sequences_by_gene_dir, g_id)
            ssbio.utils.make_dir(alignment_dir)
            protein_pickle.pairwise_align_sequences_to_representative(gapopen=gapopen, gapextend=gapextend,
                                                                      engine=engine, outdir=alignment_dir,
                                                                      parse=parse, force_rerun=force_rerun)
            protein_pickle.save_pickle(outfile=protein_seqs_aln_pickle_path)

        return g_id, protein_seqs_aln_pickle_path