Haplotype Allele Block Viewer ( python matplotlib )

NCI-H209.copynumber.tab

chr1	61735	7035990	CopyNumberRegion_0	2	0
chr1	7036801	7974076	CopyNumberRegion_1	3	0
chr1	7975788	29271013	CopyNumberRegion_2	2	0
chr1	29275883	30265652	CopyNumberRegion_3	4	2
chr1	30271899	31334415	CopyNumberRegion_4	6	2
chr1	31336903	41403766	CopyNumberRegion_5	7	2
chr1	41405663	42601214	CopyNumberRegion_6	4	2
chr1	42602991	44503629	CopyNumberRegion_7	6	2
chr1	44507772	48809180	CopyNumberRegion_8	2	0
chr1	48810111	48824773	CopyNumberRegion_9	7	1
chr1	48826008	51503719	CopyNumberRegion_10	6	2
chr1	51508246	54705383	CopyNumberRegion_11	2	0
chr1	54707046	54714379	CopyNumberRegion_12	6	0
chr1	54714820	72763371	CopyNumberRegion_13	4	2
chr1	72768418	72811149	CopyNumberRegion_14	1	0
chr1	72811903	73000937	CopyNumberRegion_15	5	0
chr1	73003082	112690327	CopyNumberRegion_16	4	2
chr1	112695673	112697391	CopyNumberRegion_17	0	0
chr1	112705956	112706150	CopyNumberRegion_18	7	0
chr1	112706194	152555528	CopyNumberRegion_19	4	2
chr1	152555706	152586541	CopyNumberRegion_20	0	0
chr1	152586576	152586577	CopyNumberRegion_21	5	0
chr1	152586594	189384373	CopyNumberRegion_22	4	2
chr1	189384741	189736833	CopyNumberRegion_23	5	2
chr1	189736933	196105810	CopyNumberRegion_24	4	2
chr1	196107264	196714874	CopyNumberRegion_25	5	2
chr1	196727850	196771281	CopyNumberRegion_26	1	0
chr1	196783626	197163536	CopyNumberRegion_27	5	0
chr1	197163927	249224389	CopyNumberRegion_28	4	2
chr2	12784	25079956	CopyNumberRegion_29	4	2
chr2	25083309	25940874	CopyNumberRegion_30	2	0
chr2	25942752	26074246	CopyNumberRegion_31	5	2
chr2	26075212	89131068	CopyNumberRegion_32	4	2
chr2	89133124	89500474	CopyNumberRegion_33	6	2
chr2	89501164	125893973	CopyNumberRegion_34	4	2
chr2	125894017	132312325	CopyNumberRegion_35	3	2
chr2	132312704	219500046	CopyNumberRegion_36	4	2
chr2	219500525	225869567	CopyNumberRegion_37	5	2
chr2	225871673	227289084	CopyNumberRegion_38	6	2
chr2	227289693	235872239	CopyNumberRegion_39	5	2
chr2	235875472	243089457	CopyNumberRegion_40	4	2
chr3	60345	53023212	CopyNumberRegion_41	2	0
chr3	53028374	53039218	CopyNumberRegion_42	0	0
chr3	53039263	90502863	CopyNumberRegion_43	2	0
chr3	93519478	93523358	CopyNumberRegion_44	8	0
chr3	93526657	97549248	CopyNumberRegion_45	6	2
chr3	97550848	98937785	CopyNumberRegion_46	5	2
chr3	98944458	98949422	CopyNumberRegion_47	0	0
chr3	98950963	102478637	CopyNumberRegion_48	5	2
chr3	102479222	107620144	CopyNumberRegion_49	6	2
chr3	107620658	114334830	CopyNumberRegion_50	7	2
chr3	114337685	130932296	CopyNumberRegion_51	4	2
chr3	130933467	131418966	CopyNumberRegion_52	2	0
chr3	131424999	131425774	CopyNumberRegion_53	5	2
chr3	131429455	174831488	CopyNumberRegion_54	4	2
chr3	174833518	176067530	CopyNumberRegion_55	2	0
chr3	176069493	176155517	CopyNumberRegion_56	5	2
chr3	176155583	197896119	CopyNumberRegion_57	4	2
chr4	12281	9859976	CopyNumberRegion_58	4	2
chr4	9864227	10214160	CopyNumberRegion_59	3	0
chr4	10223416	10235269	CopyNumberRegion_60	0	0
chr4	10236141	49658613	CopyNumberRegion_61	4	2
chr4	52685699	69372004	CopyNumberRegion_62	3	2
chr4	69374941	69489336	CopyNumberRegion_63	0	0
chr4	69505115	69506335	CopyNumberRegion_64	7	0
chr4	69521409	82339516	CopyNumberRegion_65	3	1
chr4	82347162	82386998	CopyNumberRegion_66	0	0
chr4	82390464	82413744	CopyNumberRegion_67	2	0
chr4	82415368	112146338	CopyNumberRegion_68	3	1
chr4	112149601	114359847	CopyNumberRegion_69	6	1
chr4	114359882	117367260	CopyNumberRegion_70	5	1
chr4	117371634	119165798	CopyNumberRegion_71	4	0
chr4	119166068	121855213	CopyNumberRegion_72	5	1
chr4	121855570	131577351	CopyNumberRegion_73	4	0
chr4	131581872	133871352	CopyNumberRegion_74	6	2
chr4	133873700	134005339	CopyNumberRegion_75	8	1
chr4	134009940	191027924	CopyNumberRegion_76	4	2
chr5	15532	19038761	CopyNumberRegion_77	5	2
chr5	19038937	29439275	CopyNumberRegion_78	6	2
chr5	29440909	31393022	CopyNumberRegion_79	5	2
chr5	31393075	36171420	CopyNumberRegion_80	6	3
chr5	36171932	36402794	CopyNumberRegion_81	8	3
chr5	36402904	44279496	CopyNumberRegion_82	5	2
chr5	44279545	49432832	CopyNumberRegion_83	6	2
chr5	49441958	89115932	CopyNumberRegion_84	2	0
chr5	89122730	94768511	CopyNumberRegion_85	0	0
chr5	94769498	127055318	CopyNumberRegion_86	2	0
chr5	127061322	127391993	CopyNumberRegion_87	7	0
chr5	127397909	180790321	CopyNumberRegion_88	4	0
chr6	149661	1825105	CopyNumberRegion_89	5	2
chr6	1825233	3593956	CopyNumberRegion_90	6	3
chr6	3594168	8001278	CopyNumberRegion_91	5	3
chr6	8001316	11090061	CopyNumberRegion_92	6	2
chr6	11097543	12607720	CopyNumberRegion_93	5	2
chr6	12608585	13149492	CopyNumberRegion_94	6	2
chr6	13149784	13294976	CopyNumberRegion_95	5	2
chr6	13297019	22263249	CopyNumberRegion_96	6	2
chr6	22264264	24060381	CopyNumberRegion_97	4	0
chr6	24062320	29365614	CopyNumberRegion_98	6	2
chr6	29367276	44281900	CopyNumberRegion_99	5	2
chr6	44282798	46711318	CopyNumberRegion_100	6	2
chr6	46712223	47281630	CopyNumberRegion_101	5	0
chr6	47286786	51536925	CopyNumberRegion_102	6	2
chr6	51536991	51829707	CopyNumberRegion_103	5	0
chr6	51829889	63883284	CopyNumberRegion_104	6	2
chr6	63886010	64987273	CopyNumberRegion_105	5	0
chr6	64987707	73948252	CopyNumberRegion_106	6	2
chr6	73948430	75184522	CopyNumberRegion_107	5	0
chr6	75184611	91812965	CopyNumberRegion_108	6	2
chr6	91815917	92694927	CopyNumberRegion_109	4	2
chr6	92699558	112828741	CopyNumberRegion_110	3	2
chr6	112829129	114270928	CopyNumberRegion_111	1	0
chr6	114275167	114284735	CopyNumberRegion_112	4	0
chr6	114291581	171051006	CopyNumberRegion_113	3	1
chr7	90245853	32973633	CopyNumberRegion_114	5	2
chr7	32977721	33549392	CopyNumberRegion_115	6	0
chr7	33549536	38370792	CopyNumberRegion_116	5	2
chr7	38371710	70292330	CopyNumberRegion_117	4	2
chr7	70294295	75934640	CopyNumberRegion_118	3	1
chr7	75935043	75937049	CopyNumberRegion_119	5	0
chr7	75937592	98039846	CopyNumberRegion_120	4	2
chr7	98040800	102102164	CopyNumberRegion_121	3	2
chr7	102235712	107828810	CopyNumberRegion_122	4	2
chr7	107829357	120425380	CopyNumberRegion_123	5	2
chr7	120431001	121383754	CopyNumberRegion_124	2	0
chr7	121383827	121447831	CopyNumberRegion_125	4	2
chr7	121449927	123402096	CopyNumberRegion_126	6	2
chr7	123404359	159127005	CopyNumberRegion_127	5	2
chr8	31254	3428006	CopyNumberRegion_128	7	2
chr8	3429068	3589281	CopyNumberRegion_129	6	0
chr8	3589415	14790841	CopyNumberRegion_130	7	2
chr8	14793512	15410250	CopyNumberRegion_131	6	2
chr8	15410393	15709310	CopyNumberRegion_132	5	0
chr8	15709998	18501993	CopyNumberRegion_133	6	2
chr8	18502306	18866338	CopyNumberRegion_134	5	0
chr8	18866367	19056987	CopyNumberRegion_135	6	2
chr8	19060538	24974477	CopyNumberRegion_136	4	2
chr8	24974522	24984333	CopyNumberRegion_137	0	0
chr8	24991104	24992008	CopyNumberRegion_138	5	2
chr8	24993669	61689295	CopyNumberRegion_139	4	2
chr8	61692829	61731073	CopyNumberRegion_140	8	1
chr8	61732263	68330527	CopyNumberRegion_141	4	2
chr8	68336207	70127366	CopyNumberRegion_142	6	2
chr8	70132518	71145689	CopyNumberRegion_143	5	2
chr8	71146049	118660792	CopyNumberRegion_144	6	2
chr8	118671466	119819044	CopyNumberRegion_145	5	2
chr8	119819845	123779599	CopyNumberRegion_146	6	2
chr8	123780071	132208969	CopyNumberRegion_147	5	2
chr8	132209925	133798940	CopyNumberRegion_148	6	2
chr8	133799279	137071329	CopyNumberRegion_149	5	2
chr8	137073679	142164446	CopyNumberRegion_150	6	2
chr8	142165340	146298156	CopyNumberRegion_151	5	2
chr9	40909	141091395	CopyNumberRegion_152	4	2
chr10	72759	17246240	CopyNumberRegion_153	3	1
chr10	17246584	17254832	CopyNumberRegion_154	5	2
chr10	17255588	19544126	CopyNumberRegion_155	4	2
chr10	19545383	42433753	CopyNumberRegion_156	3	2
chr10	42433815	42469720	CopyNumberRegion_157	5	2
chr10	42469893	135506705	CopyNumberRegion_158	4	2
chr11	198509	4968117	CopyNumberRegion_159	4	2
chr11	4968510	4983161	CopyNumberRegion_160	7	0
chr11	4985908	92083681	CopyNumberRegion_161	4	2
chr11	92089406	92640674	CopyNumberRegion_162	6	2
chr11	92640786	134944770	CopyNumberRegion_163	4	2
chr12	150442	889010	CopyNumberRegion_164	4	2
chr12	889638	1306656	CopyNumberRegion_165	2	0
chr12	1306953	1350676	CopyNumberRegion_166	5	2
chr12	1354344	31280591	CopyNumberRegion_167	4	2
chr12	31280903	31409879	CopyNumberRegion_168	6	2
chr12	31410283	133778190	CopyNumberRegion_169	4	2
chr13	19026949	115108398	CopyNumberRegion_170	4	0
chr14	19002124	106224525	CopyNumberRegion_171	4	2
chr14	106246302	106863527	CopyNumberRegion_172	6	2
chr14	106864722	107285437	CopyNumberRegion_173	4	2
chr15	20016328	102469041	CopyNumberRegion_174	4	2
chr16	60777	2178576	CopyNumberRegion_175	4	2
chr16	2232794	3645909	CopyNumberRegion_176	2	0
chr16	3649400	3878722	CopyNumberRegion_177	0	0
chr16	3887250	3887251	CopyNumberRegion_178	5	0
chr16	3888921	67951233	CopyNumberRegion_179	4	2
chr16	67952968	70779144	CopyNumberRegion_180	2	0
chr16	70785196	71202798	CopyNumberRegion_181	3	1
chr16	71205194	72424159	CopyNumberRegion_182	2	0
chr16	72428052	72868230	CopyNumberRegion_183	5	2
chr16	72868464	90287536	CopyNumberRegion_184	4	2
chr17	526	60381847	CopyNumberRegion_185	4	0
chr17	60382144	61252362	CopyNumberRegion_186	5	0
chr17	61257548	80572628	CopyNumberRegion_187	4	0
chr17	80574531	81048659	CopyNumberRegion_188	3	0
chr18	11542	1134703	CopyNumberRegion_189	9	2
chr18	1142463	3800808	CopyNumberRegion_190	6	3
chr18	3807348	30335775	CopyNumberRegion_191	5	2
chr18	30336716	32240117	CopyNumberRegion_192	6	2
chr18	32241543	68218793	CopyNumberRegion_193	5	2
chr18	68222284	71197081	CopyNumberRegion_194	6	2
chr18	71198887	78015057	CopyNumberRegion_195	5	2
chr19	90910	59097855	CopyNumberRegion_196	4	2
chr20	61305	54903467	CopyNumberRegion_197	4	2
chr20	54904078	55643296	CopyNumberRegion_198	3	1
chr20	55644151	55645708	CopyNumberRegion_199	6	0
chr20	55646251	62054969	CopyNumberRegion_200	4	2
chr20	62056012	62956154	CopyNumberRegion_201	3	1
chr21	10736871	48096958	CopyNumberRegion_202	4	2
chr22	16052528	25661376	CopyNumberRegion_203	4	2
chr22	25661697	25918709	CopyNumberRegion_204	2	0
chr22	25922342	25922348	CopyNumberRegion_205	7	2
chr22	25922548	51234456	CopyNumberRegion_206	4	2
chrX	168477	2693624	CopyNumberRegion_207	6	2
chrX	2694239	7023493	CopyNumberRegion_208	4	0
chrX	7024973	88465611	CopyNumberRegion_209	2	0
chrX	88465982	92421465	CopyNumberRegion_210	3	1
chrX	92421505	154929499	CopyNumberRegion_211	2	0
chrX	154943318	155254182	CopyNumberRegion_212	4	1
chrY	179541	58977523	CopyNumberRegion_213	2	0
chrY	58977827	59017393	CopyNumberRegion_214	3	0
chrY	59017524	59018406	CopyNumberRegion_215	0	0
chrY	59019879	59303517	CopyNumberRegion_216	3	0

pHapBlockBinomialTest.step1.sh

 for i in H3K27ac H3K36me3 H3K4me1 H3K4me3 H3K9me3
  do
 #     for j in chr1 chr2 chr3 chr4 chr5 chr6 chr7 chr8 chr9 chr10 chr11 chr12 chr13 chr14 chr15 chr16 chr17 chr18 chr19 chr20 chr21 chr22 chrX chrY

#     do
# ./xHapBlockBinomialTest.py ./HAP/test_$j.haplotype ../Mapping/${i}_H209_hg19_mapping_v2.sorted.bam ./NCI-H209.copynumber.tab> ./HapBlockBinomialTest/${i}.HapBlockBinomialTest.tab.$j
#  grep Simple ./HapBlockBinomialTest/$i.HapBlockBinomialTest.tab.$j > ./HapBlockBinomialTest/Simple/$i.HapBlockBinomialTest.tab.$j.simple
#  
#
# done
 cat ./HapBlockBinomialTest/Simple/$i.*.simple > ./HapBlockBinomialTest/Raw/$i.HapBlockBinomialTest.raw.tab
 done

xMergeHapBlockTestToTable.step2.py

#!/usr/bin/python
# programmer : zhuxp
# usage:
class HapSimple:
    def __init__(self,x):
        self.chr=x[0]
        self.start=int(x[1])# 0 index
        self.stop=int(x[2])
        self.id=x[3]
        self.coverage=int(x[6])
        self.phase0=int(x[7])
        self.phase1=int(x[8])
        self.alleleA=int(x[9])
        self.alleleB=int(x[10])
        self.p1=float(x[11])
        self.p2=float(x[12])
    def __cmp__(self,other):
        return cmp(self.chr,other.chr) or cmp(self.start,other.start)
def Main():
    hHap=dict()
    lName=["H3K27ac","H3K36me3","H3K4me1","H3K4me3","H3K9me3"]
    fName=[]
    print "chr\tstart\tstop\tid\tcopynumbers\t","\t".join(lName)
    for i in lName:
        fName.append(i+".HapBlockBinomialTest.raw.tab")
    for i,fname in enumerate(fName):
        f=open(fname)
        for line in f:
            line=line.strip()
            if line[0]=="#": continue
            if len(line)==0: continue
            x=line.split("\t")
            a=HapSimple(x)
            if hHap.has_key(a.id):
                hHap[a.id][lName[i]]=a
            else:
                hHap[a.id]=dict()
                hHap[a.id][lName[i]]=a
    ids=hHap.keys()
    ids.sort(key=lambda item:hHap[item]["H3K27ac"] )

    for id in ids:
        a=hHap[id]
        print a["H3K27ac"].chr,"\t",a["H3K27ac"].start,"\t",a["H3K27ac"].stop,"\t",id,"\t",
        print str(a["H3K27ac"].alleleA)+":"+str(a["H3K27ac"].alleleB),
        flag=0
        for i in lName:
            print "\t",str(a[i].phase0)+":"+str(a[i].phase1),
            if (a[i].p1<0.1 and a[i].p2<0.1):
                if flag==0:
                    print "*",
                    flag=cmp(a[i].phase0,a[i].phase1)
                else:
                    if flag==cmp(a[i].phase0,a[i].phase1):
                        print "*",
                    else:
                        print "-*",

                #if(a[i].phase0 >= a[i].phase1): print "+",
                #else : print "-",
        print


if __name__=="__main__":
    Main()

xViewHapblock.step3.py

#!/usr/bin/python
# programmer : zhuxp
# usage: sys.argv[0]  haploblock.tab  gene.tab chr1:1-100000 
'''
 draw figures of haploblock and genes
 input file : gene.tab (download from UCSC table)
              haploblock distribution file ( Generated from ReadRaw.py and *.simple table)
                            *.simple table were generated from xHaploBlockBinomialTest.py (pHaploBlockBinomialTest.sh)
                            the input for that program is BAM file and HapCUT result file
'''

'''
KNOWN ISSUES: No gridspec in sysbio
'''

import sys
import numpy as np
import matplotlib.pyplot as plt

import matplotlib.gridspec as gridspec
from zlab.zbed import *

def addGeneToFig(gene,ax,name=0):
    '''
    add gene to figures
    from green to red is the gene direction
    '''
    if name==1:
        ax.text(gene.start,15,gene.id)
    cds=gene.new_cds()
    utr5=gene.new_utr5()
    utr3=gene.new_utr3()
    cds_exons=cds.Exons()
    for cds_exon in cds_exons:
        ax.bar(cds_exon.start,10,cds_exon.stop-cds_exon.start,color="blue",alpha=0.5)
        ax.bar(cds_exon.start,-10,cds_exon.stop-cds_exon.start,color="blue",alpha=0.5)
    if not utr3 is None:
        for utr3_exon in utr3.Exons():
            ax.bar(utr3_exon.start,5,utr3_exon.stop-utr3_exon.start,color="red",alpha=0.5)
            ax.bar(utr3_exon.start,-5,utr3_exon.stop-utr3_exon.start,color="red",alpha=0.5)
    if not utr5 is None:
        for utr5_exon in utr5.Exons():
            ax.bar(utr5_exon.start,5,utr5_exon.stop-utr5_exon.start,color="green",alpha=0.5)
            ax.bar(utr5_exon.start,-5,utr5_exon.stop-utr5_exon.start,color="green",alpha=0.5)
    for intron in gene.Introns():
        ax.bar(intron.start,1,intron.stop-intron.start,color="yellow",alpha=0.5)
        ax.bar(intron.start,-1,intron.stop-intron.start,color="yellow",alpha=0.5)
        #print intron


def parse(x):
    '''
    parse x:y in Blocks
    '''
    a=x.split(" ")
    b=a[0].split(":")
    c=[0,0,"."]
    c[0]=int(b[0])
    c[1]=int(b[1])
    if len(a)>1:
        c[2]=a[1]
    return c

class tBlock(Bed):
    def __init__(self,x):
        self.chr=x[0].rstrip()
        self.start=int(x[1])
        self.stop=int(x[2])
        self.id=x[3]
        self.copy_numbers=parse(x[4])
        self.hEpi=dict()
        self.hEpi["H3K27ac"]=parse(x[5])
        self.hEpi["H3K36me3"]=parse(x[6])
        self.hEpi["H3K4me1"]=parse(x[7])
        self.hEpi["H3K4me3"]=parse(x[8])
        self.hEpi["H3K9me3"]=parse(x[9])
        self.allele_flag=0
        for i in self.hEpi.keys():
            a=self.hEpi[i]
            if a[2]=="*":
                self.allele_flag=1




def addToFig(block,ax,pos,j,d):
        k=block.hEpi.keys()
        k.sort()
        phase=[list(),list()]
        for i in k:
            a=block.hEpi[i]
            phase[0].append(a[0])
            phase[1].append(-a[1])
        print phase[0]
        print phase[1]
        print k
        ind=np.arange(5)
        ind+=pos
        width=1
        rects0=ax.bar(ind,phase[0],width,color=["r","g","b","y","black"])
        rects1=ax.bar(ind,phase[1],width,color=["r","g","b","y","black"])
        ax.text((j+0.5)/21,0.8,"%d"%d,transform=ax.transAxes) 
        for i,x in enumerate(k):
            a=block.hEpi[x]
            if a[2]=="*" and a[0]>a[1]:
                 height = rects0[i].get_height()
                 ax.text(rects0[i].get_x()+rects0[i].get_width()/2., 1.05*height, '*',ha='center', va='bottom')
            if a[2]=="*" and a[1]>a[0]:
                 height = rects1[i].get_height()
                 ax.text(rects1[i].get_x()+rects0[i].get_width()/2., -1.05*height, '*',ha='center', va='top')
        

        
            
                
        
        
def region_parse(r):
    r=r.strip()
    a=r.split(":")
    b=a[1].split("-")
    return (a[0],int(b[0]),int(b[1]))

def Main():
    region=sys.argv[3]
    (region_chr,region_start,region_stop)=region_parse(sys.argv[3])
    regionBed=Bed([region_chr,region_start-1,region_stop,"QueryRegion",".",0])
    '''
    READING GENE into LIST
    '''
    genes=list()
    genefile=open(sys.argv[2])
    for line in genefile:
        line=line.strip()
        if line[0]=="#": continue
        if len(line)==0: continue
        gene=GeneBed(line.split("\t"))
        if Bed.overlap(gene,regionBed):
            genes.append(gene)
            print gene
    genes.sort()
    '''
    READING BLOCK INTO LIST
    '''
    f=open(sys.argv[1])
    blocks=list()
    for line in f:
        line=line.strip()
        if line[0]=="#": continue
        if len(line)==0: continue
        x=line.split("\t")
        if not x[1].isdigit: continue
        try:
            
            block=tBlock(x)
        except:
            continue
        if Bed.overlap(block,regionBed):
            blocks.append(block)
    blocks.sort()

    
    '''
    initialize fig and draw the first panel of chromosome 
    '''
    number=len(blocks)
    panel_number=number/20+1
    if number%20==0: panel_number-1
    fig = plt.figure(figsize=(25,panel_number*4+3), frameon = True )
    ratios=[4 for i in range(panel_number+2)]
    ratios[0]=2
    ratios[1]=1
    gs = gridspec.GridSpec(panel_number+2, 1, width_ratios=[1],height_ratios=ratios)
    '''
    Draw Blocks Distribution in Chromosome in the first panel
    '''
    ax1=plt.subplot(gs[1])
    ax1.set_yticks([])
    ax1.axis([regionBed.start,regionBed.stop,0,10])
    ax1.set_ylim(0,10)
    ax1.set_xticks([])
    for i,x in enumerate(blocks):
        if x.allele_flag==1:
            ax1.bar(x.start,10,x.stop-x.start,color="r")
        else:
            ax1.bar(x.start,10,x.stop-x.start,color="b")
        if (i+1)%5==0: 
            height=-3
            ax1.text((x.start+x.stop)/2,height,'%d'%(i+1),ha='center',va='bottom')
    '''
    Draw Genes Distribution in the second panel
    '''
    ax2=plt.subplot(gs[0])
    ax2.text(0.05,0.7,region_chr+":"+str(region_start)+"-"+str(region_stop),transform=ax2.transAxes,color="black") 
    ax2.set_yticks([])
    ax2.axis([regionBed.start,regionBed.stop,-30,30])
    for i in genes:
        addGeneToFig(i,ax2,1)
    
    '''
    Draw blocks in other panels 
    '''
    #gsBlocks=gs[2:]
    ax=[0 for i in range(panel_number)]

    for i in range(panel_number):
        ax[i]=plt.subplot(gs[i+2])
    for i,x in enumerate(blocks):
        addToFig(x,ax[i/20],(i%20)*10+5,i%20,i+1)
        ax[i/20].set_xticks([])    
        ax[i/20].set_xlim(0,210)
    
    ax[-1].text(0.95,0.9,"H3K27ac",transform=ax[-1].transAxes,color="r") 
    ax[-1].text(0.95,0.7,"H3K36me3",transform=ax[-1].transAxes,color="g") 
    ax[-1].text(0.95,0.5,"H3K4me1",transform=ax[-1].transAxes,color="b") 
    ax[-1].text(0.95,0.3,"H3K4me3",transform=ax[-1].transAxes,color="y") 
    ax[-1].text(0.95,0.1,"H3K9me3",transform=ax[-1].transAxes,color="black") 
 #   ax[-1].text(0.95,0.9,["H3K27ac","H3K36me3","H3K4me1","H3K4me3","H3K9me3"],transform=ax[-1].transAxes,color=["r","g","b","y","black"]) 
 #   ax[-1].legend()
    '''
    Show Figures
    '''
    plt.show()    
    fig.savefig(region_chr+"_"+str(region_start)+"_"+str(region_stop)+".png")


    
if __name__=="__main__":
    Main()