Property which dynamically construct transcript objects for all
        transcript IDs associated with this gene.
6916
def transcripts(self):
        
        transcript_id_results = self.db.query(
            select_column_names=['transcript_id'],
            filter_column='gene_id',
            filter_value=self.id,
            feature='transcript',
            distinct=False,
            required=False)

        # We're doing a SQL query for each transcript ID to fetch
        # its particular information, might be more efficient if we
        # just get all the columns here, but how do we keep that modular?
        return [
            self.genome.transcript_by_id(result[0])
            for result in transcript_id_results
        ]
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
249872
def from_gene_ids_and_names(cls, gene_names: Dict[str, str]):
        """"""
        genes = [ExpGene(id_, name=name) for id_, name in gene_names.items()]
        return cls.from_genes(genes)
1235044
def _collapse_transcripts(in_file, window, data, out_dir, include_gene_names=True):
    
    if out_dir is None:
        out_dir = os.path.dirname(in_file)
    out_file = os.path.join(out_dir,
                            "%s-transcripts_w%s.bed" % (os.path.splitext(os.path.basename(in_file))[0],
                                                        window))
    chrom_sizes = {}
    for contig in ref.file_contigs(dd.get_ref_file(data), data["config"]):
        chrom_sizes[contig.name] = contig.size
    if not utils.file_uptodate(out_file, in_file):
        with file_transaction(data, out_file) as tx_out_file:
            prep_file = "%s-sortprep%s" % os.path.splitext(tx_out_file)
            sort_cmd = bedutils.get_sort_cmd(os.path.dirname(tx_out_file))
            cmd = "{sort_cmd} -k4,4 -k1,1 {in_file} > {prep_file}"
            do.run(cmd.format(**locals()), "Sort BED file by transcript name")
            with open(tx_out_file, "w") as out_handle:
                # Work around for segmentation fault issue with groupby
                # https://github.com/daler/pybedtools/issues/131#issuecomment-89832476
                x = pybedtools.BedTool(prep_file)
                def gen():
                    for r in x:
                        yield r
                for name, rs in itertools.groupby(gen(), lambda r: (r.name, r.chrom)):
                    rs = list(rs)
                    r = rs[0]
                    if r.chrom in chrom_sizes:
                        for gcoords in _group_coords(rs):
                            min_pos = max(min(gcoords) - window, 0)
                            max_pos = min(max(gcoords) + window, chrom_sizes[r.chrom])
                            if include_gene_names:
                                out_handle.write("%s\t%s\t%s\t%s\n" % (r.chrom, min_pos, max_pos, r.name))
                            else:
                                out_handle.write("%s\t%s\t%s\n" % (r.chrom, min_pos, max_pos))
    return bedutils.sort_merge(out_file, data)
391371
def exon_by_id(self, exon_id):
        
        if exon_id not in self._exons:
            field_names = [
                "seqname",
                "start",
                "end",
                "strand",
                "gene_name",
                "gene_id",
            ]

            contig, start, end, strand, gene_name, gene_id = self.db.query_one(
                select_column_names=field_names,
                filter_column="exon_id",
                filter_value=exon_id,
                feature="exon",
                distinct=True)

            self._exons[exon_id] = Exon(
                exon_id=exon_id,
                contig=contig,
                start=start,
                end=end,
                strand=strand,
                gene_name=gene_name,
                gene_id=gene_id)

        return self._exons[exon_id]
254664
def get_genomic_seq_for_transcript(self, transcript_id, expand):
        
        
        headers = {"content-type": "application/json"}
        
        self.attempt = 0
        ext = "/sequence/id/{0}?type=genomic;expand_3prime={1};expand_5prime={1}".format(transcript_id, expand)
        r = self.ensembl_request(ext, headers)
        
        gene = json.loads(r)
        
        seq = gene["seq"]
        seq_id = gene["id"]
        
        if seq_id != transcript_id:
            raise ValueError("ensembl gave the wrong transcript")
        
        desc = gene["desc"].split(":")
        chrom = desc[2]
        start = int(desc[3]) + expand
        end = int(desc[4]) - expand
        strand_temp = int(desc[5])
        
        strand = "+"
        if strand_temp == -1:
            strand = "-"
        
        return (chrom, start, end, strand, seq)
766764
def get_transcript_preferences(course_id):
    
    try:
        transcript_preference = TranscriptPreference.objects.get(course_id=course_id)
    except TranscriptPreference.DoesNotExist:
        return

    return TranscriptPreferenceSerializer(transcript_preference).data
1000213
def cdsparse(self, record):
        
        try:
            # Find genes that are present in all strains of interest - the number of times the gene is found is
            # equal to the number of strains. Earlier parsing ensures that the same gene is not present in a strain
            # more than once
            if self.genes[self.genenames[record.id]] == len(self.runmetadata.samples):
                # Add the gene names and sequences to the appropriate dictionaries
                try:
                    self.genesequence[self.genenames[record.id]].add(str(record.seq))
                # Initialise the dictionary as required, then populate as above
                except KeyError:
                    self.genesequence[self.genenames[record.id]] = set()
                    self.genesequence[self.genenames[record.id]].add(str(record.seq))
                try:
                    self.coresequence[str(record.seq)].add(record.id)
                except KeyError:
                    self.coresequence[str(record.seq)] = set()
                    self.coresequence[str(record.seq)].add(record.id)
        except KeyError:
            pass
1365235
def transcript_id_of_protein_id(self, protein_id):
        
        results = self._query_transcript_ids(
            "protein_id",
            protein_id,
            feature="CDS")
        if len(results) == 0:
            raise ValueError("Protein ID not found: %s" % protein_id)
        assert len(results) == 1, \
            ("Should have only one transcript ID for a given protein ID, "
             "but found %d: %s" % (len(results), results))
        return results[0]
1127586
def gtf(args):
    
    p = OptionParser(gtf.__doc__)
    opts, args = p.parse_args(args)

    if len(args) != 1:
        sys.exit(not p.print_help())

    gffile, = args
    gff = Gff(gffile)
    transcript_info = AutoVivification()
    for g in gff:
        if g.type.endswith(("RNA", "transcript")):
            if "ID" in g.attributes and "Parent" in g.attributes:
                transcript_id = g.get_attr("ID")
                gene_id = g.get_attr("Parent")
            elif "mRNA" in g.attributes and "Gene" in g.attributes:
                transcript_id = g.get_attr("mRNA")
                gene_id = g.get_attr("Gene")
            else:
                transcript_id = g.get_attr("ID")
                gene_id = transcript_id
            transcript_info[transcript_id]["gene_id"] = gene_id
            transcript_info[transcript_id]["gene_type"] = g.type
            continue

        if g.type not in valid_gff_to_gtf_type.keys():
            continue

        try:
            transcript_id = g.get_attr("Parent", first=False)
        except IndexError:
            transcript_id = g.get_attr("mRNA", first=False)

        g.type = valid_gff_to_gtf_type[g.type]
        for tid in transcript_id:
            if tid not in transcript_info: continue
            gene_type = transcript_info[tid]["gene_type"]
            if not gene_type.endswith("RNA") and not gene_type.endswith("transcript"):
                continue
            gene_id = transcript_info[tid]["gene_id"]
            g.attributes = dict(gene_id=[gene_id], transcript_id=[tid])
            g.update_attributes(gtf=True, urlquote=False)

            print(g)
613990
def ids_from_seqs_iterative(seqs, app, query_parser, \
    scorer=keep_everything_scorer, max_iterations=None, blast_db=None,\
    max_seqs=None, ):
    
    if isinstance(scorer, int):
        scorer = make_shotgun_scorer(scorer)
    seqs_to_check = list(seqs)
    checked_ids = {}
    curr_iteration = 0
    while seqs_to_check:
        unchecked_ids = {}
        #pass seqs to command
        all_output = app(seqs_to_check)
        output = all_output.get('BlastOut', all_output['StdOut'])

        for query_id, match_id, match_score in query_parser(output):
            if query_id not in checked_ids:
                checked_ids[query_id] = {}
            checked_ids[query_id][match_id] = match_score
            if match_id not in checked_ids:
                unchecked_ids[match_id] = True
        all_output.cleanUp()
        if unchecked_ids:
            seq_file = fasta_cmd_get_seqs(unchecked_ids.keys(),
                app.Parameters['-d'].Value)['StdOut']
            seqs_to_check = []
            for s in FastaCmdFinder(fasta_cmd_get_seqs(\
                unchecked_ids.keys(), app.Parameters['-d'].Value)['StdOut']):
                seqs_to_check.extend(s)
        else:
            seqs_to_check = []
        #bail out if max iterations or max seqs was defined and we've reached it
        curr_iteration += 1
        if max_iterations and (curr_iteration >= max_iterations):
            break
        if max_seqs:
            curr = scorer(checked_ids)
            if len(curr) >= max_seqs:
                return curr
    return scorer(checked_ids)
860718
def _get_genes(self, variant):
        
        ensembl_ids = []
        hgnc_symbols = []
        
        for transcript in variant.transcripts:
            if transcript.ensembl_id:
                ensembl_ids.append(transcript.ensembl_id)
            if transcript.hgnc_symbol:
                hgnc_symbols.append(transcript.hgnc_symbol)
        
        genes = get_gene_info(
                        ensembl_ids=ensembl_ids, 
                        hgnc_symbols=hgnc_symbols
                        )
        return genes
305381
def generate_transcript_counts(data):
    """"""
    data["count_file"] = featureCounts.count(data)

    if dd.get_fusion_mode(data, False) and not dd.get_fusion_caller(data):
        oncofuse_file = oncofuse.run(data)
        if oncofuse_file:
            data = dd.set_oncofuse_file(data, oncofuse_file)

    if dd.get_transcriptome_align(data):
        # to create a disambiguated transcriptome file realign with bowtie2
        if dd.get_disambiguate(data):
            logger.info("Aligning to the transcriptome with bowtie2 using the "
                        "disambiguated reads.")
            bam_path = data["work_bam"]
            fastq_paths = alignprep._bgzip_from_bam(bam_path, data["dirs"], data, is_retry=False, output_infix='-transcriptome')
            if len(fastq_paths) == 2:
                file1, file2 = fastq_paths
            else:
                file1, file2 = fastq_paths[0], None
            ref_file = dd.get_ref_file(data)
            data = bowtie2.align_transcriptome(file1, file2, ref_file, data)
        else:
            file1, file2 = dd.get_input_sequence_files(data)
        if not dd.get_transcriptome_bam(data):
            ref_file = dd.get_ref_file(data)
            logger.info("Transcriptome alignment was flagged to run, but the "
                        "transcriptome BAM file was not found. Aligning to the "
                        "transcriptome with bowtie2.")
            data = bowtie2.align_transcriptome(file1, file2, ref_file, data)
    data = spikein.counts_spikein(data)
    return [[data]]
80608
def cli(self, prt=sys.stdout):
        """"""
        kws = self.objdoc.get_docargs(prt=None)
        if os.path.exists(kws['i']):
            obj = NCBIgeneFileReader(kws['i'])
            nts = obj.get_nts()
            if nts:
                geneid2nt = self._get_geneid2nt(nts)
                self._wrpy_ncbi_gene_nts(kws['o'], geneid2nt, prt)
        else:
            raise RuntimeError("\n{DOC}\n**ERROR: NO FILE FOUND: {NCBI}".format(
                NCBI=kws['i'], DOC=__doc__))
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();

		}

	}
347685
def canonical_transcripts(gtf, out_file):
    
    if file_exists(out_file):
        return out_file
    db = get_gtf_db(gtf)
    with file_transaction(out_file) as tx_out_file:
        with open(tx_out_file, "w") as out_handle:
            for gene in db.features_of_type('gene'):
                exon_list = []
                for ti, transcript in enumerate(db.children(gene, level=1)):
                    cds_len = 0
                    total_len = 0
                    exons = list(db.children(transcript, level=1))
                    for exon in exons:
                        exon_length = len(exon)
                        if exon.featuretype == 'CDS':
                            cds_len += exon_length
                        total_len += exon_length

                    exon_list.append((cds_len, total_len, transcript, exons))

                # If we have CDS, then use the longest coding transcript
                if max(i[0] for i in exon_list) > 0:
                    best = sorted(exon_list)[0]
                # Otherwise, just choose the longest
                else:
                    best = sorted(exon_list, key=lambda x: x[1])[0]
                for exon in db.children(best[2], level=1):
                    out_handle.write(str(exon) + "\n")
    return out_file
142964
def index_transcriptome(gtf_file, ref_file, data):
    
    gtf_fasta = gtf.gtf_to_fasta(gtf_file, ref_file)
    bowtie2_index = os.path.splitext(gtf_fasta)[0]
    bowtie2_build = config_utils.get_program("bowtie2", data["config"]) + "-build"
    cmd = "{bowtie2_build} --offrate 1 {gtf_fasta} {bowtie2_index}".format(**locals())
    message = "Creating transcriptome index of %s with bowtie2." % (gtf_fasta)
    do.run(cmd, message)
    return bowtie2_index