Python GenomicRegionSet.gene_association方法代码示例

# 需要导入模块: from rgt.GenomicRegionSet import GenomicRegionSet [as 别名]
# 或者: from rgt.GenomicRegionSet.GenomicRegionSet import gene_association [as 别名]
def get_dbss(input_BED,output_BED,rna_fasta,output_rbss,organism,l,e,c,fr,fm,of,mf,rm,temp):
    regions = GenomicRegionSet("Target")
    regions.read_bed(input_BED)
    regions.gene_association(organism=organism, show_dis=True)

    connect_rna(rna_fasta, temp=temp, rna_name="RNA")
    rnas = SequenceSet(name="rna", seq_type=SequenceType.RNA)
    rnas.read_fasta(os.path.join(temp,"rna_temp.fa"))
    rna_regions = get_rna_region_str(os.path.join(temp,rna_fasta))
    # print(rna_regions)
    genome = GenomeData(organism)
    genome_path = genome.get_genome()
    txp = find_triplex(rna_fasta=rna_fasta, dna_region=regions, 
                       temp=temp, organism=organism, remove_temp=False,
                       l=l, e=e, c=c, fr=fr, fm=fm, of=of, mf=mf, genome_path=genome_path,
                       prefix="targeted_region", dna_fine_posi=True)

    print("Total binding events:\t",str(len(txp)))
    txp.write_bed(output_BED)
    txp.write_txp(filename=output_BED.replace(".bed",".txp"))
    rbss = txp.get_rbs()
    dbd_regions(exons=rna_regions, sig_region=rbss, rna_name="rna", output=output_rbss, 
                out_file=True, temp=temp, fasta=False)

# 需要导入模块: from rgt.GenomicRegionSet import GenomicRegionSet [as 别名]
# 或者: from rgt.GenomicRegionSet.GenomicRegionSet import gene_association [as 别名]
def rna_associated_gene(rna_regions, name, organism):
    if rna_regions:
        s = [ rna_regions[0][0], min([e[1] for e in rna_regions]), 
              max([e[2] for e in rna_regions]), rna_regions[0][3] ]
        g = GenomicRegionSet("RNA associated genes")
        g.add( GenomicRegion(chrom=s[0], initial=s[1], final=s[2], name=name, orientation=s[3]) )
        asso_genes = g.gene_association(organism=organism, promoterLength=1000, show_dis=True)

        genes = asso_genes[0].name.split(":")
        closest_genes = []
        for n in genes:
            if name not in n: closest_genes.append(n)
        closest_genes = set(closest_genes)

        if len(closest_genes) == 0:
            return "."
        else:
            return ":".join(closest_genes)
    else:
        return "."

# 需要导入模块: from rgt.GenomicRegionSet import GenomicRegionSet [as 别名]
# 或者: from rgt.GenomicRegionSet.GenomicRegionSet import gene_association [as 别名]
        bed1 = GenomicRegionSet("input")
        bed1.read_bed(args.i)
        bed2 = bed1.mergebyname()
        bed2.write_bed(args.o)

    ############### BED rename regions #######################################
    elif args.mode == "bed_rename":
        print(tag+": [BED] Rename regions by associated genes")
        print("input:\t" + args.i)
        print("output:\t" + args.o)
        print("organism:\t" + args.organism)

        bed = GenomicRegionSet(args.i)
        bed.read_bed(args.i)
        renamebed = bed.gene_association(gene_set=None, organism=args.organism, 
                                         promoterLength=args.l, 
                                         threshDist=args.t, show_dis=args.d)
        renamebed.write_bed(args.o)

    ############### BED extend ###############################################
    elif args.mode == "bed_extend":
        print("input:\t" + args.i)
        print("output:\t" + args.o)
        bed = GenomicRegionSet("bed")
        bed.read_bed(args.i)
        for region in bed:
            if args.oz:
                if region.initial == region.final:
                    region.final += args.l
            else:
                if args.both:

# 需要导入模块: from rgt.GenomicRegionSet import GenomicRegionSet [as 别名]
# 或者: from rgt.GenomicRegionSet.GenomicRegionSet import gene_association [as 别名]
class RandomTest:
    def __init__(self, rna_fasta, rna_name, dna_region, organism, showdbs=False):
        self.organism = organism
        genome = GenomeData(organism)
        self.genome_path = genome.get_genome()
        # RNA: Path to the FASTA file
        self.rna_fasta = rna_fasta
        self.showdbs = showdbs

        rnas = SequenceSet(name="rna", seq_type=SequenceType.RNA)
        rnas.read_fasta(self.rna_fasta)
        if rna_name:
            self.rna_name = rna_name
        else:
            self.rna_name = rnas[0].name

        # DNA: GenomicRegionSet
        self.dna_region = GenomicRegionSet(name="target")
        self.dna_region.read_bed(dna_region)
        self.dna_region = self.dna_region.gene_association(organism=self.organism, show_dis=True)

        self.topDBD = []
        self.stat = OrderedDict(name=rna_name, genome=organism)
        self.stat["target_regions"] = str(len(self.dna_region))


    def get_rna_region_str(self, rna):
        """Getting the rna region from the information header with the pattern:
                REGION_chr3_51978050_51983935_-_"""
        self.rna_regions = get_rna_region_str(rna)
        if self.rna_regions and len(self.rna_regions[0]) == 5:
            self.rna_expression = float(self.rna_regions[0][-1])
        else:
            self.rna_expression = "n.a."

    def connect_rna(self, rna, temp):
        d = connect_rna(rna, temp, self.rna_name)
        self.stat["exons"] = str(d[0])
        self.stat["seq_length"] = str(d[1])
        self.rna_len = d[1]

    def target_dna(self, temp, remove_temp, cutoff, l, e, c, fr, fm, of, mf, par, obed=False):
        """Calculate the true counts of triplexes on the given dna regions"""
        self.triplexator_p = [ l, e, c, fr, fm, of, mf ]

        txp = find_triplex(rna_fasta=os.path.join(temp, "rna_temp.fa"), dna_region=self.dna_region,
                           temp=temp, organism=self.organism, remove_temp=remove_temp,
                           l=l, e=e, c=c, fr=fr, fm=fm, of=of, mf=mf, par=par, genome_path=self.genome_path,
                           prefix="targeted_region", dna_fine_posi=False)
        txp.merge_rbs(rm_duplicate=True, region_set=self.dna_region, asgene_organism=self.organism, cutoff=cutoff)
        self.txp = txp
        self.stat["DBSs_target_all"] = str(len(self.txp))
        txp.remove_duplicates()
        self.rbss = txp.merged_dict.keys()
        # if len(self.rbss) == 0:
        #     print("ERROR: No potential binding event. Please change the parameters.")
        #     sys.exit(1)

        txpf = find_triplex(rna_fasta=os.path.join(temp, "rna_temp.fa"), dna_region=self.dna_region,
                            temp=temp, organism=self.organism, remove_temp=remove_temp,
                            l=l, e=e, c=c, fr=fr, fm=fm, of=of, mf=mf, par=par, genome_path=self.genome_path,
                            prefix="dbs", dna_fine_posi=True)
        txpf.remove_duplicates()
        txpf.merge_rbs(rbss=self.rbss, rm_duplicate=True, asgene_organism=self.organism)
        self.txpf = txpf

        self.stat["DBSs_target_all"] = str(len(self.txpf))

        self.counts_tr = OrderedDict()
        self.counts_dbs = OrderedDict()

        for rbs in self.rbss:
            tr = len(self.txp.merged_dict[rbs])
            self.counts_tr[rbs] = [tr, len(self.dna_region) - tr]
            self.counts_dbs[rbs] = len(self.txpf.merged_dict[rbs])

        self.region_dbd = self.txpf.sort_rbs_by_regions(self.dna_region)

        self.region_dbs = self.txpf.sort_rd_by_regions(regionset=self.dna_region)
        self.region_dbsm = {}
        self.region_coverage = {}

        for region in self.dna_region:
            self.region_dbsm[region.toString()] = self.region_dbs[region.toString()].get_dbs().merge(w_return=True)
            self.region_coverage[region.toString()] = float(self.region_dbsm[region.toString()].total_coverage()) / len \
                (region)
        self.stat["target_regions"] = str(len(self.dna_region))

        if obed:
            # btr = self.txp.get_dbs()
            # btr = btr.gene_association(organism=self.organism, show_dis=True)
            # btr.write_bed(os.path.join(temp, obed + "_target_region_dbs.bed"))
            # dbss = txpf.get_dbs()
            # dbss.write_bed(os.path.join(temp, obed + "_dbss.bed"))

            # output = self.dna_region.gene_association(organism=self.organism, show_dis=True)

            self.txp.write_bed(filename=os.path.join(temp, obed + "_target_region_dbs.bed"),
                               dbd_tag=False,
                               remove_duplicates=False, associated=self.organism)
#.........这里部分代码省略.........