Mercurial > repos > big-tiandm > sirna_plant
diff chr_plot.pl @ 0:07745c0958dd draft
Uploaded
| author | big-tiandm |
|---|---|
| date | Thu, 18 Sep 2014 21:40:25 -0400 |
| parents | |
| children |
line wrap: on
line diff
--- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/chr_plot.pl Thu Sep 18 21:40:25 2014 -0400 @@ -0,0 +1,607 @@ +#!/usr/bin/perl -w +#========================================================================================== +# Date: +# Title: +# Comment: Program to plot gene structure +# Input: 1. input file of Gene region annotation which format like GenePred +# 2. input file of Transcripts region annotation which format like GenePred +# 3. input file of gene snp detail info +# Output: output file of gene structure graph by html or svg formt +# Test Usage: +#======================================================================================== +#use strict; +my $version=1.00; +use SVG; +use Getopt::Long; +my %opt; +GetOptions(\%opt,"g=s","l=s","chro=s","mark=s","span=s","te=s","t=s","cen=s","c=s","o=s","out=s","h"); +if (!(defined $opt{g} and defined $opt{c} and defined $opt{l} and defined $opt{chro} and defined $opt{mark} and defined $opt{o}) || defined $opt{h}) { +&usage; +} +my $span=$opt{span}; +#===============================Define Attribute========================================== +my %attribute=( + canvas=>{ + 'width'=>1500, + 'height'=>1800 + }, + text=>{ + 'stroke'=>"#000000", + 'fill'=>"none", + 'stroke-width'=>0.5 + }, + line=>{ + 'stroke'=>"black", + 'stroke-width'=>1 + }, + csv=>{ + 'stroke'=>"red", + 'stroke-width'=>0.5 + }, + exon=>{ + 'stroke'=>"black", + 'stroke-width'=>1 + }, + intron=>{ + 'stroke'=>"black", + 'stroke-width'=>1.5 + }, + font=>{ + 'fill'=>"#000000", + 'font-size'=>12, + 'font-size2'=>10, + #'font-weight'=>'bold', + 'font-family'=>"Arial" + #'font-family'=>"ArialNarrow-bold" + }, + rect=>{ + 'fill'=>"lightgreen", + 'stroke'=>"black", + 'stroke-width'=>0.5 + }, + readwidth=>0.5 +); +#############################s#define start coordinate and scale + +my $length=$opt{"l"};#11 +#my @target_e_value=qw(1.78 1.83 1.92 2.00 1.92 2.00 1.93 2.11 2.05 2.03 1.92 1.91 1.54 1.72 1.67); + +my $chr=$opt{"chro"}; +my $start=1; +my $XOFFSET=50; +my $YOFFSET=60; +#my $length=$end-$start+1; +my $Xscale=600/$length;#定义X轴比例尺 1:1000 x轴的坐标长度都要按照此比例尺换算 +#my $high_cov=$high_cov9B1=0.5;#定义峰图最高峰 +#my $Yscale=1/$high_cov;#定义Y轴比例尺 1:60 y轴的坐标长度都要按照此比例尺换算 +#========================================New canvas============================ +#--------------------------------------------------------------- +if (defined($opt{cen})) { + open(CEN,"$opt{cen}")||die"cannot open the file $opt{cen}"; + my %centromere; + while (my $aline=<CEN>) { + chomp $aline; + next if($aline=~/^\#/); + my @temp=split/\t/,$aline; + $centromere{$temp[0]}[0]=$temp[1]; + $centromere{$temp[0]}[1]=$temp[2]; + } + close CEN; +} + +#### Starting #### +#新建画布 +my $svg=SVG->new(); +#画图起始点 +my $canvas_start_x=$XOFFSET; +my $canvas_end_x=$XOFFSET+$length*$Xscale;#按照比例尺 画线 +my $canvas_start_y=$YOFFSET; +my $canvas_end_y=$YOFFSET; +#Draw a straight line between two points start(x1,y1) and end(x2,y2). +#location attribute +$svg->line(id=>'l1',x1=>$canvas_start_x,y1=>$canvas_start_y,x2=>$canvas_end_x,y2=>$canvas_end_y,'stroke',$attribute{line}{'stroke'},'stroke-width',$attribute{line}{'stroke-width'}); +$long_scale=int ($length/10);#十等分 大刻度 +#大坐标刻度 +for ($i=0;$i<=10;$i++) { + my $long_x_start=$XOFFSET+$long_scale*$i*$Xscale; + my $long_x_end=$long_x_start; + my $long_y_start=$YOFFSET; + my $long_y_end=$YOFFSET-5; + $svg->line('x1',$long_x_start,'y1',$long_y_start,'x2',$long_x_end,'y2',$long_y_end,'stroke',$attribute{line}{'stroke'},'stroke-width',$attribute{line}{'stroke-width'}); + my $Bscale=$start+$long_scale*$i; + my $cdata=int ($Bscale/1000000); + $svg->text('x',$long_x_start,'y',$long_y_start-10,'style','fill:black;text-anchor:middle','stroke',$attribute{text}{'stroke'},'stroke-width',$attribute{text}{'stroke-width'},'font-size',12,'font-family',$attribute{font}{'font-family'},'-cdata',$cdata."M"); +} + +if (defined($opt{cen})) { + $svg->rect('x',$XOFFSET+$centromere{$chr2}[0]*$Xscale,'y',$YOFFSET-2,'width',($centromere{$chr2}[1]-$centromere{$chr2}[0]+1)*$Xscale,'height',5,'stroke',"black",'stroke-width',"5",'fill',"black"); + + $svg->rect('x',$XOFFSET+$length*$Xscale+15,'y',$YOFFSET+20,'width',10,'height',5,'stroke',"black",'stroke-width',$attribute{intron}{'stroke-width'},'fill',"black"); + $svg->text('x',$XOFFSET+$length*$Xscale+35,'y',$YOFFSET+20+5,'stroke',"black",'stroke-width',$attribute{text}{'stroke-width'},'font-size',$attribute{font}{'font-size'},'font-family',$attribute{font}{'font-family'},'-cdata',"Centromere"); +} + +$svg->text('x',$XOFFSET+$length*$Xscale*0.42,'y',$YOFFSET-30,'stroke',"black",'stroke-width',1,'font-size',15,'font-family',$attribute{font}{'font-family'},'-cdata',$chr); +#====================================================MAIN PROGRAM======================== + +#===========================================gene list data================================ +#open (TXT,">$opt{out}"); +#open(TE,"$opt{te}")||die"cannot open the file $opt{te}"; +#my %te; +#while (my $aline=<TE>) { +# chomp $aline; +# my @temp=split/\t/,$aline; +# if ($temp[2] eq "Y") { +# $te{$temp[0]}=1; +# } +#} +#close TE; +#===================================Exp gene list data ================================= +#open(TARGET,"$opt{t}")||die"cannot open the file $opt{t}"; +#my %target_e; +#my %target_rpkm; +#while (my $aline=<TARGET>) { +# chomp $aline;##Gene 19B1 1PA1 1LC1 29B2 2PA2 2LC2 39B3 3PA3 3LC3 93-4C PA-4C LY-4C 9343 PA43 LY43 19B1_VS_1LC1 19B1_VS_1PA1 1PA1_VS_1LC1 29B2_VS_2LC2 29B2_VS_2PA2 2PA2_VS_2LC2 39B3_VS_3LC3 39B3_VS_3PA3 3PA3_VS_3LC3 934C_LY4C 934C_PA4C PA4C_LY4C 9343_VS_LY43 9343_VS_PA43 PA43_VS_LY43 mpv_1 fold_1 mpv_2 fold_2 mpv_3 fold_3 mpv_4_B fold_4_B mpv_leaf fold_lead +# next if($aline=~/^\#/); +# my @temp=split/\t/,$aline; +# $target_rpkm{$temp[0]}="$temp[7]\t$temp[8]\t$temp[9]"; +# if ($temp[7]>$target_e_value[6]||$temp[8]>$target_e_value[7]||$temp[9]>$target_e_value[8]) { +# $target_e{$temp[0]}="$temp[7]\t$temp[8]\t$temp[9]"; +# } +#} +#close TARGET; +#==================================================================================== +my @genelist; +#my @te_genelist; +my @target_e; +open(GENE,"$opt{g}")||die"cannot open the file $opt{g}"; +while (my $aline=<GENE>) { + chomp $aline;#LOC_Os01g01280 Chr1 133291 134685 + + my @temp=split/\t/,$aline; + #my $te; + if ($temp[1]=~/^Chr(\d)$/) { + $temp[1]="Chr0$1"; + } + next if($temp[1] ne $chr); + #push @genelist,[$temp[0],$temp[3],$temp[4]]; + push @genelist,[$temp[0],$temp[2],$temp[3]]; +# if ($te{$temp[0]}) { +# push @te_genelist,[$temp[0],$temp[3],$temp[4]]; +# } +# else{ +# push @genelist,[$temp[0],$temp[3],$temp[4]]; +# } +# if ($target_e{$temp[0]}) { +# push @target_e,[$temp[0],$temp[3],$temp[4]]; +# } + +} +close GENE; +my @gene_desity; +#my @region_target_rpkm; +@genelist=sort{$a->[1] <=> $b->[1]}@genelist; +for (my $i=0;$i<@genelist ;$i++) { +# if ($genelist[$i][1]>($num1+1)*$span) { +# $num1++; +# } +# if ($genelist[$i][1]>$num1*$span&&$genelist[$i][1]<=($num1+1)*$span) { +# $gene_desity[$num1]++; +# +# } + my $start=int($genelist[$i][1]/$span); + my $end=int($genelist[$i][2]/$span); + #my @t_rpkm=split/\t/,$target_rpkm{$genelist[$i][0]}; + if ($start==$end) { + $gene_desity[$start]++; + #for (my $rs=0;$rs<3 ;$rs++) { + #print TXT "$rs\t$i\t$genelist[$i][0]\t$target_rpkm{$genelist[$i][0]}\n"; + #$region_target_rpkm[$start][$rs]+=$t_rpkm[$rs]; + #} + + #$target_rpkm_desity[$start][0]+=$temp[0]; + #$target_rpkm_desity[$start][1]+=$temp[1]; + #$target_rpkm_desity[$start][2]+=$temp[2]; + } + else{ + for (my $k=$start;$k<=$end ;$k++) { + $gene_desity[$k]++; + #for (my $rs=0;$rs<3 ;$rs++) { + #$region_target_rpkm[$k][$rs]+=$t_rpkm[$rs]; + #} + #$target_rpkm_desity[$k][0]+=$temp[0]; + #$target_rpkm_desity[$k][1]+=$temp[1]; + #$target_rpkm_desity[$k][2]+=$temp[2]; + } + } +} +#------------------------------------------region_gene_number------------------------- +my $max_gene_number=0; +my $total=0; +for (my $i=0;$i<@gene_desity ;$i++) { + if (!(defined($gene_desity[$i]))) { + $gene_desity[$i]=0; + } + if ($gene_desity[$i]>$max_gene_number) { + $max_gene_number=$gene_desity[$i]; + } + #print TXT "$i\t$gene_desity[$i]\n"; + $total+=$gene_desity[$i]; +} +print "Gene max:$max_gene_number\ntotal:$total\n"; +my $abc=@genelist; +#my $cba=@te_genelist; +print "aaaaaa:$abc\n"; +#print "bbbbbb:$cba\n"; +#---------------------------------------------region_gene_rpkm------------------------ +#my $max_region_gene_rpkm=0; +#my @max_region_gene_rpkm=qw(0 0 0); +#my @total_region_gene_rpkm=qw(0 0 0); +#for (my $i=0;$i<@region_target_rpkm ;$i++) { +# for (my $s=0; $s<3;$s++) { +# +# if (!(defined($region_target_rpkm[$i][$s]))) { +# $region_target_rpkm[$i][$s]=0; +# } +# $total_region_gene_rpkm[$s]+=$region_target_rpkm[$i][$s]; +# if ($max_region_gene_rpkm<$region_target_rpkm[$i][$s]) { +# $max_region_gene_rpkm=$region_target_rpkm[$i][$s]; +# } +# if ($max_region_gene_rpkm[$s]<$region_target_rpkm[$i][$s]) { +# $max_region_gene_rpkm[$s]=$region_target_rpkm[$i][$s]; +# } +# } +#} +#my @ave_region_gene_rpkm=qw(0 0 0); +#my $max_ave_rpkm=0; +#for (my $i=0;$i<3;$i++) { +# $ave_region_gene_rpkm[$i]=$total_region_gene_rpkm[$i]/($#region_target_rpkm+1); +# if ($max_ave_rpkm<$ave_region_gene_rpkm[$i]) { +# $max_ave_rpkm=$ave_region_gene_rpkm[$i]; +# } +#} +# +#print "***max region gene rpkm :$max_region_gene_rpkm\n\n"; +#print "@max_region_gene_rpkm\n"; +#if ($max_region_gene_rpkm>10*$max_ave_rpkm) { +# $max_region_gene_rpkm=10*$max_ave_rpkm; +#} +#my @max_region_rpkm; +#---------------------------------------------------------------------------------- +#my @te_gene_desity; +#@te_genelist=sort{$a->[1] <=> $b->[1]}@te_genelist; +##my $num2=0; +#for (my $i=0;$i<@te_genelist ;$i++) { +## if ($te_genelist[$i][1]>($num2+1)*$span) { +## $num2=int($te_genelist[$i][1]/$span); +## } +## if ($te_genelist[$i][1]>$num2*$span&&$te_genelist[$i][1]<=($num2+1)*$span) { +## $te_gene_desity[$num2]++; +## } +# my $start=int($te_genelist[$i][1]/$span); +# my $end=int($te_genelist[$i][2]/$span); +# if ($start==$end) { +# $te_gene_desity[$start]++; +# } +# else{ +# for (my $k=$start;$k<=$end ;$k++) { +# $te_gene_desity[$k]++; +# } +# } +#} +#$max_te_gene_number=0; +#$total=0; +#for (my $i=0;$i<@te_gene_desity ;$i++) { +# if (!(defined($te_gene_desity[$i]))) { +# $te_gene_desity[$i]=0; +# } +# if ($te_gene_desity[$i]>$max_te_gene_number) { +# $max_te_gene_number=$te_gene_desity[$i]; +# } +# print TXT "$i\t$te_gene_desity[$i]\n"; +# $total+=$te_gene_desity[$i]; +#} +# +#print "TE gene max:$max_te_gene_number\ntotal:$total\n"; +#------------------------------------------------------- +#my @target_e_desity; +#@target_e=sort{$a->[1] <=> $b->[1]}@target_e; +#for (my $i=0;$i<@target_e ;$i++) { +# my $start=int($target_e[$i][1]/$span); +# my $end=int($target_e[$i][2]/$span); +# if ($start==$end) { +# $target_e_desity[$start]++; +# } +# else{ +# for (my $k=$start;$k<=$end ;$k++) { +# $target_e_desity[$k]++; +# } +# } +#} +#my $max_target_e_number=0; +#$total=0; +#for (my $i=0;$i<@target_e_desity ;$i++) { +# if (!(defined($target_e_desity[$i]))) { +# $target_e_desity[$i]=0; +# } +# if ($target_e_desity[$i]>$max_target_e_number) { +# $max_target_e_number=$target_e_desity[$i]; +# } +# print TXT "$i\t$target_e_desity[$i]\n"; +# $total+=$target_e_desity[$i]; +#} +# +#print "Target max:$max_target_e_number\ntotal:$total\n"; + +#====================================cluster data======================================= +open(CLUSTER,"$opt{c}")||die"cannot open the file $opt{c}"; +my $mark="$opt{mark}"; +my @sample=split/\#/,$mark; +my $sample_num=@sample; +my @cluster; +my @cluster_density; +my @max_cluster_read=(0)x$sample_num; +while (my $aline=<CLUSTER>) { + next if($aline=~/^\#/); + chomp $aline;##ID chr strand start end 19B1 + #clusterID major_length percent sample1 sample2 sample3 + #Chr01:192429-192452 24nt 1.00 0.00 97222.22 0.00 + my @temp=split/\t/,$aline; + my @id=split/\:/,$temp[0]; + my @posit=split/\-/,$id[1]; + next if($id[0] ne $chr);#Chr.01 chr04 + push @cluster,[$id[0],$posit[0],$posit[1],@temp[3..3+$sample_num+1]]; + #push @cluster,[$temp[0],$temp[3],$temp[4],$temp[11],$temp[12],$temp[13]];#ID start end rpkm(19B1,1PA1,1LC1); + for (my $i=0;$i<@sample ;$i++) { + if($temp[3+$i]>$max_cluster_read[$i]){ + $max_cluster_read[$i]=$temp[3+$i]; + } + } +} +close CLUSTER; +@cluster=sort{$a->[1] <=> $b->[1]}@cluster; +for(my $i=0;$i<$#cluster;$i++) { + my $start=int($cluster[$i][1]/$span); + my $end=int($cluster[$i][2]/$span); + if ($start==$end) { + for (my $j=0;$j<$sample_num ;$j++) { + $cluster_density[$start][$j]+=$cluster[$i][$j+$sample_num]; + } + + } + else{ + for (my $m=$start;$m<=$end ;$m++) { + for (my $j=0;$j<$sample_num ;$j++) { + $cluster_density[$m][$j]+=$cluster[$i][$j+$sample_num]; + } + } + } +} +my @max_cluster_density=(0)x$sample_num; +my @total_cluster_density=(0)x$sample_num; +my $max_cluster_density=0; +for (my $i=0;$i<@sample;$i++) { + for (my $k=0;$k<$#cluster_density ;$k++) { + + if (!(defined($cluster_density[$k][$i]))) { + $cluster_density[$k][$i]=0; + } + $total_cluster_density[$i]+=$cluster_density[$k][$i]; + if ($cluster_density[$k][$i]>$max_cluster_density[$i]) { + $max_cluster_density[$i]=$cluster_density[$k][$i]; + } + if ($cluster_density[$k][$i]>$max_cluster_density) { + $max_cluster_density=$cluster_density[$k][$i]; + + } + } +} +my @ave_cluster_density=(0)x$sample_num; +my $max_ave=0; +for (my $i=0;$i<@sample;$i++) { + $ave_cluster_density[$i]=$total_cluster_density[$i]/($#cluster_density+1); + if ($max_ave<$ave_cluster_density[$i]) { + $max_ave=$ave_cluster_density[$i]; + } +} + +print "max cluster reads:@max_cluster_read\n"; +print "max cluster density:@max_cluster_density\n"; +if ($max_cluster_density>10*$max_ave) { + $max_cluster_density=10*$max_ave; +} +#===================================== plot gene desity ================================= +#===row +my $Y1scale=5; +my $y1_gene_density=$YOFFSET+$max_gene_number*$Y1scale+10; +for (my $i=0;$i<@gene_desity-1 ;$i++) { + my $density_start_x=$XOFFSET+$i*$span*$Xscale; + my $density_start_y=$y1_gene_density-$gene_desity[$i]*$Y1scale; + my $density_end_x=$XOFFSET+($i+1)*$span*$Xscale; + my $density_end_y=$y1_gene_density-$gene_desity[$i+1]*$Y1scale; + my $heigth=$gene_desity[$i]/$max_gene_number; + $svg->rect('x',$density_start_x,'y',$density_start_y,'width',$span*$Xscale,'height',$y1_gene_density-$density_start_y,'stroke',"read",'stroke-width',$attribute{intron}{'stroke-width'},'fill',"green"); + #$svg->line('x1',$density_start_x,'y1',$density_start_y,'x2',$density_end_x,'y2',$density_end_y,'stroke',"red",'stroke-width',$attribute{csv}{'stroke-width'}); + +} +$svg->line('x1',$XOFFSET,'y1',$y1_gene_density,'x2',$XOFFSET+$length*$Xscale,'y2',$y1_gene_density,'stroke',"black",'stroke-width',"black"); + +#====clomun +$svg->line('x1',$XOFFSET,'y1',$y1_gene_density,'x2',$XOFFSET,'y2',$y1_gene_density-$max_gene_number*$Y1scale,'stroke',"black",'stroke-width',"black"); +$svg->text('x',$XOFFSET-15,'y',$y1_gene_density,'stroke',"black",'stroke-width',$attribute{text}{'stroke-width'},'font-size',$attribute{font}{'font-size2'},'font-family',$attribute{font}{'font-family'},'-cdata',"0");#min gene number +$svg->text('x',$XOFFSET-20,'y',$y1_gene_density-$max_gene_number*$Y1scale,'stroke',"black",'stroke-width',$attribute{text}{'stroke-width'},'font-size',$attribute{font}{'font-size2'},'font-family',$attribute{font}{'font-family'},'-cdata',$max_gene_number);#max gene number +#========================================================================================= +#=============================plot TE gene desity========================================= +#===row +=cut +my $Y2scale=$Y1scale; +my $y2_te_gene_density=$y1_gene_density; +for (my $i=0;$i<@te_gene_desity-1 ;$i++) { + my $te_density_start_x=$XOFFSET+$i*$span*$Xscale; + my $te_density_start_y=$y2_te_gene_density+$te_gene_desity[$i]*$Y2scale; + my $te_density_end_x=$XOFFSET+($i+1)*$span*$Xscale; + my $te_density_end_y=$y2_te_gene_density+$te_gene_desity[$i+1]*$Y2scale; + #my $te_heigth=$te_gene_desity[$i]/$max_gene_number; + $svg->rect('x',$te_density_start_x,'y',$y2_te_gene_density,'width',$span*$Xscale,'height',$te_density_start_y-$y2_te_gene_density,'stroke',"read",'stroke-width',$attribute{intron}{'stroke-width'},'fill',"green"); + #$svg->line('x1',$density_start_x,'y1',$density_start_y,'x2',$density_end_x,'y2',$density_end_y,'stroke',"red",'stroke-width',$attribute{csv}{'stroke-width'}); + +} +#column +$svg->line('x1',$XOFFSET,'y1',$y2_te_gene_density,'x2',$XOFFSET,'y2',$y2_te_gene_density+$max_te_gene_number*$Y2scale,'stroke',"black",'stroke-width',"black"); +$svg->text('x',$XOFFSET-20,'y',$y2_te_gene_density+$max_te_gene_number*$Y2scale,'stroke',"black",'stroke-width',$attribute{text}{'stroke-width'},'font-size',$attribute{font}{'font-size2'},'font-family',$attribute{font}{'font-family'},'-cdata',$max_te_gene_number);#min gene number +=cut +#=======================gene density txt================================================== +my $md=$span/1000; +#$svg->text('x',$XOFFSET-30,'y',$YOFFSET+10,'width',"698",'height',"298",'stroke',"black",'stroke-width',$attribute{text}{'stroke-width'},'font-size',$attribute{font}{'font-size2'},'font-family',$attribute{font}{'font-family'},'-cdata',"Number per \\30 kb","rotate","90",'writing-mode',"tb-rl"); + +$svg->rect('x',$XOFFSET+$length*$Xscale+15,'y',$y1_gene_density-$max_gene_number*$Y1scale*0.7,'width',10,'height',5,'stroke',"green",'stroke-width',$attribute{intron}{'stroke-width'},'fill',"green"); +$svg->text('x',$XOFFSET+$length*$Xscale+35,'y',$y1_gene_density-$max_gene_number*$Y1scale*0.7+5,'stroke',"green",'stroke-width',$attribute{text}{'stroke-width'},'font-size',$attribute{font}{'font-size'},'font-family',$attribute{font}{'font-family'},'-cdata',"Genes"); + +#$svg->rect('x',$XOFFSET+$length*$Xscale+15,'y',$y1_gene_density-$max_gene_number*$Y1scale*0.7+20,'width',10,'height',5,'stroke',"green",'stroke-width',$attribute{intron}{'stroke-width'},'fill',"green"); +#$svg->text('x',$XOFFSET+$length*$Xscale+35,'y',$y1_gene_density-$max_gene_number*$Y1scale*0.7+20+5,'stroke',"green",'stroke-width',$attribute{text}{'stroke-width'},'font-size',$attribute{font}{'font-size'},'font-family',$attribute{font}{'font-family'},'-cdata',"TE Genes"); + +$svg->text('x',$XOFFSET+$length*$Xscale+35,'y',$y1_gene_density-$max_gene_number*$Y1scale*0.7+20,'stroke',"black",'stroke-width',$attribute{text}{'stroke-width'},'font-size',$attribute{font}{'font-size'},'font-family',$attribute{font}{'font-family'},'-cdata',"genes number \/ $md kb"); +#========================================================================================= +#=============================plot exp gene desity========================================= +=cut +my $y3_target_e_density=$y2_te_gene_density+$max_te_gene_number*$Y2scale+$max_target_e_number*$Y2scale+20; +#my $Y3scale=$Y1scale; +for (my $i=0;$i<@target_e_desity-1 ;$i++) { + my $target_e_density_start_x=$XOFFSET+$i*$span*$Xscale; + my $target_e_density_start_y=$y3_target_e_density-$target_e_desity[$i]*$Y2scale; + my $target_e_density_end_x=$XOFFSET+($i+1)*$span*$Xscale; + my $target_e_density_end_y=$y3_target_e_density-$target_e_desity[$i+1]*$Y2scale; + $svg->rect('x',$target_e_density_start_x,'y',$target_e_density_start_y,'width',$span*$Xscale,'height',$y3_target_e_density-$target_e_density_start_y,'stroke',"read",'stroke-width',$attribute{intron}{'stroke-width'},'fill',"red"); + +} +$svg->line('x1',$XOFFSET,'y1',$y3_target_e_density,'x2',$XOFFSET+$length*$Xscale,'y2',$y3_target_e_density,'stroke',"black",'stroke-width',"black"); +#column +$svg->line('x1',$XOFFSET,'y1',$y3_target_e_density,'x2',$XOFFSET,'y2',$y3_target_e_density-$max_te_gene_number*$Y2scale,'stroke',"black",'stroke-width',"black"); +$svg->text('x',$XOFFSET-15,'y',$y3_target_e_density,'stroke',"black",'stroke-width',$attribute{text}{'stroke-width'},'font-size',$attribute{font}{'font-size2'},'font-family',$attribute{font}{'font-family'},'-cdata',0);#max gene number +$svg->text('x',$XOFFSET-20,'y',$y3_target_e_density-$max_te_gene_number*$Y2scale+10,'stroke',"black",'stroke-width',$attribute{text}{'stroke-width'},'font-size',$attribute{font}{'font-size2'},'font-family',$attribute{font}{'font-family'},'-cdata',$max_target_e_number);#max gene number +#=========================================exp gene indensity txt========================== +$svg->rect('x',$XOFFSET+$length*$Xscale+15,'y',$y3_target_e_density-$max_target_e_number*$Y2scale*0.7,'width',10,'height',5,'stroke',"red",'stroke-width',$attribute{intron}{'stroke-width'},'fill',"red"); +$svg->text('x',$XOFFSET+$length*$Xscale+35,'y',$y3_target_e_density-$max_target_e_number*$Y2scale*0.7+5,'stroke',"red",'stroke-width',$attribute{text}{'stroke-width'},'font-size',$attribute{font}{'font-size'},'font-family',$attribute{font}{'font-family'},'-cdata',"Exp Genes"); +$svg->text('x',$XOFFSET+$length*$Xscale+35,'y',$y3_target_e_density-$max_target_e_number*$Y2scale*0.7+20,'stroke',"black",'stroke-width',$attribute{text}{'stroke-width'},'font-size',$attribute{font}{'font-size'},'font-family',$attribute{font}{'font-family'},'-cdata',"genes number \/ 50kb"); +#calculate the different +=cut +#======================================================================================== +my $Y3scale=0.04/$max_cluster_density*2500; +#my $y3_cluster_density=$y2_te_gene_density+$max_te_gene_number*5+$max_cluster_density*$Y3scale+10; +my $y3_cluster_density=$y1_gene_density+20; +my $y3_total=$y1_gene_density+10; +my @cluster_color=qw(fuchsia violet tomato); +for (my $s=0;$s<3 ;$s++) { + $y3_total=$y3_total+$max_cluster_density*$Y3scale+5; + $svg->line('x1',$XOFFSET,'y1',$y3_total,'x2',$XOFFSET+$length*$Xscale,'y2',$y3_total,'stroke',"black",'stroke-width',$attribute{csv}{'stroke-width'}); + $svg->line('x1',$XOFFSET,'y1',$y3_total,'x2',$XOFFSET,'y2',$y3_total-$max_cluster_density*$Y3scale,'stroke',"black",'stroke-width',$attribute{csv}{'stroke-width'}); + + $svg->text('x',$XOFFSET-15,'y',$y3_total,'stroke',"black",'stroke-width',$attribute{text}{'stroke-width'},'font-size',$attribute{font}{'font-size2'},'font-family',$attribute{font}{'font-family'},'-cdata',"0");#min gene number + + if ($max_cluster_density[$s]>$max_cluster_density) { + $svg->text('x',$XOFFSET-40,'y',$y3_total-$max_cluster_density*$Y3scale+10,'stroke',"black",'stroke-width',$attribute{text}{'stroke-width'},'font-size',$attribute{font}{'font-size2'},'font-family',$attribute{font}{'font-family'},'-cdata',$max_cluster_density[$s]);#max gene number + } + else{ + $svg->text('x',$XOFFSET-40,'y',$y3_total-$max_cluster_density[$s]*$Y3scale+10,'stroke',"black",'stroke-width',$attribute{text}{'stroke-width'},'font-size',$attribute{font}{'font-size2'},'font-family',$attribute{font}{'font-family'},'-cdata',$max_cluster_density[$s]);#max gene number + } + + for (my $i=0;$i<$#cluster_density ;$i++) { + my $cluster_density_start_x=$XOFFSET+$i*$span*$Xscale; + my $cluster_density_end_x=$XOFFSET+($i+1)*$span*$Xscale; + if ($cluster_density[$i][$s]>$max_cluster_density) { + $cluster_density[$i][$s]=$max_cluster_density; + } + if ($cluster_density[$i+1][$s]>$max_cluster_density) { + $cluster_density[$i+1][$s]=$max_cluster_density; + } + my $cluster_density_start_y=$y3_total-$cluster_density[$i][$s]*$Y3scale; + my $cluster_density_end_y=$y3_total-$cluster_density[$i+1][$s]*$Y3scale; + $svg->line('x1',$cluster_density_start_x,'y1',$cluster_density_start_y,'x2',$cluster_density_end_x,'y2',$cluster_density_end_y,'stroke',$cluster_color[$s],'stroke-width',$attribute{csv}{'stroke-width'}); + } +} +for (my $s=0;$s<@sample ;$s++) { + $svg->line('x1',$XOFFSET+$length*$Xscale+10,'y1',$y3_cluster_density+($y3_total-$y3_cluster_density)*0.3+30*$s,'x2',$XOFFSET+$length*$Xscale+30,'y2',$y3_cluster_density+($y3_total-$y3_cluster_density)*0.3+30*$s,'stroke',$cluster_color[$s],'stroke-width',$attribute{csv}{'stroke-width'}); + $svg->text('x',$XOFFSET+$length*$Xscale+35,'y',$y3_cluster_density+($y3_total-$y3_cluster_density)*0.3+30*$s+5,'stroke',$cluster_color[$s],'stroke-width',$attribute{text}{'stroke-width'},'font-size',$attribute{font}{'font-size'},'font-family',$attribute{font}{'font-family'},'-cdata',"Small RNAs ".$sample[$s]); +} + +$svg->text('x',$XOFFSET+$length*$Xscale+35,'y',$y3_cluster_density+($y3_total-$y3_cluster_density)*0.3+30*3-5,'stroke',"black",'stroke-width',$attribute{text}{'stroke-width'},'font-size',$attribute{font}{'font-size'},'font-family',$attribute{font}{'font-family'},'-cdata',"indensity of sRNA \/ $md kb"); +#===================================plot region target gene rpkm================ +=cut +my $y4_region_gene_rpkm_1=$y3_total+10; +my $y4_region_gene_rpkm=$y4_region_gene_rpkm_1; +my $Y4scale=0.1/$max_region_gene_rpkm*1000; +my @cluster_color_t=qw(blue slateblue steelblue); +for (my $s=0;$s<3 ;$s++) { + $y4_region_gene_rpkm+=$max_region_gene_rpkm*$Y4scale+5; + + $svg->line('x1',$XOFFSET,'y1',$y4_region_gene_rpkm,'x2',$XOFFSET+$length*$Xscale,'y2',$y4_region_gene_rpkm,'stroke',"black",'stroke-width',$attribute{csv}{'stroke-width'}); + + $svg->line('x1',$XOFFSET,'y1',$y4_region_gene_rpkm,'x2',$XOFFSET,'y2',$y4_region_gene_rpkm-$max_region_gene_rpkm*$Y4scale,'stroke',"black",'stroke-width',$attribute{csv}{'stroke-width'}); + + $svg->text('x',$XOFFSET-15,'y',$y4_region_gene_rpkm,'stroke',"black",'stroke-width',$attribute{text}{'stroke-width'},'font-size',$attribute{font}{'font-size2'},'font-family',$attribute{font}{'font-family'},'-cdata',"0");#min gene number + + if ($max_region_gene_rpkm[$s]>$max_region_gene_rpkm) { + $svg->text('x',$XOFFSET-40,'y',$y4_region_gene_rpkm-$max_region_gene_rpkm*$Y4scale+10,'stroke',"black",'stroke-width',$attribute{text}{'stroke-width'},'font-size',$attribute{font}{'font-size2'},'font-family',$attribute{font}{'font-family'},'-cdata',$max_region_gene_rpkm[$s]);#max gene number + } + else{ + $svg->text('x',$XOFFSET-40,'y',$y4_region_gene_rpkm-$max_region_gene_rpkm[$s]*$Y4scale+10,'stroke',"black",'stroke-width',$attribute{text}{'stroke-width'},'font-size',$attribute{font}{'font-size2'},'font-family',$attribute{font}{'font-family'},'-cdata',$max_region_gene_rpkm[$s]);#max gene number + } + for (my $i=0;$i<$#region_target_rpkm ;$i++) { + my $region_target_rpkm_start_x=$XOFFSET+$i*$span*$Xscale; + my $region_target_rpkm_end_x=$XOFFSET+($i+1)*$span*$Xscale; + if ($region_target_rpkm[$i][$s]>$max_region_gene_rpkm) { + $region_target_rpkm[$i][$s]=$max_region_gene_rpkm; + } + if ($region_target_rpkm[$i+1][$s]>$max_region_gene_rpkm) { + $region_target_rpkm[$i+1][$s]=$max_region_gene_rpkm; + } + my $region_target_rpkm_start_y=$y4_region_gene_rpkm-$region_target_rpkm[$i][$s]*$Y4scale; + my $region_target_rpkm_end_y=$y4_region_gene_rpkm-$region_target_rpkm[$i+1][$s]*$Y4scale; + $svg->line('x1',$region_target_rpkm_start_x,'y1',$region_target_rpkm_start_y,'x2',$region_target_rpkm_end_x,'y2',$region_target_rpkm_end_y,'stroke',$cluster_color_t[$s],'stroke-width',$attribute{csv}{'stroke-width'}); + } + +} +for (my $s=0;$s<3 ;$s++) { + $svg->line('x1',$XOFFSET+$length*$Xscale+10,'y1',$y4_region_gene_rpkm_1+($y4_region_gene_rpkm-$y4_region_gene_rpkm_1)*0.3+30*$s,'x2',$XOFFSET+$length*$Xscale+30,'y2',$y4_region_gene_rpkm_1+($y4_region_gene_rpkm-$y4_region_gene_rpkm_1)*0.3+30*$s,'stroke',$cluster_color_t[$s],'stroke-width',$attribute{csv}{'stroke-width'}); + $svg->text('x',$XOFFSET+$length*$Xscale+35,'y',$y4_region_gene_rpkm_1+($y4_region_gene_rpkm-$y4_region_gene_rpkm_1)*0.3+30*$s+5,'stroke',$cluster_color_t[$s],'stroke-width',$attribute{text}{'stroke-width'},'font-size',$attribute{font}{'font-size'},'font-family',$attribute{font}{'font-family'},'-cdata',"Target Genes ".$sample[$s]); +} +$svg->text('x',$XOFFSET+$length*$Xscale+35,'y',$y4_region_gene_rpkm_1+($y4_region_gene_rpkm-$y4_region_gene_rpkm_1)*0.3+30*3-5,'stroke',"black",'stroke-width',$attribute{text}{'stroke-width'},'font-size',$attribute{font}{'font-size'},'font-family',$attribute{font}{'font-family'},'-cdata',"indensity of genes \/ 50kb"); +=cut +#for (my $s=0;$s<3 ;$s++) { +# $svg->line('x1',$XOFFSET+$length*$Xscale+10,'y1',$y4_region_gene_rpkm_1+($y4_region_gene_rpkm-$y3_cluster_density)*0.3+30*$s,'x2',$XOFFSET+$length*$Xscale+30,'y2',$y4_region_gene_rpkm_1+($y4_region_gene_rpkm-$y3_cluster_density)*0.3+30*$s,'stroke',$cluster_color[$s],'stroke-width',$attribute{csv}{'stroke-width'}); +# $svg->text('x',$XOFFSET+$length*$Xscale+35,'y',$y4_region_gene_rpkm_1+($y4_region_gene_rpkm-$y3_cluster_density)*0.3+30*$s+5,'stroke',$cluster_color[$s],'stroke-width',$attribute{text}{'stroke-width'},'font-size',$attribute{font}{'font-size'},'font-family',$attribute{font}{'font-family'},'-cdata',$sample[$s]); +#} +#========================================================================================= + +#sub ExpG_number { +# +#} +#sub ExpCluster_number{ +# +#} + + +#======================================================================================== +open (OUT,">$opt{o}"); +print OUT $svg->xmlify(); + +sub log2 { + my $n = shift; + return log($n)/log(2); +} + +sub usage{ +print <<"USAGE"; +Version $version +Usage: +$0 +options: +-l centromere length +-chro +-mark \# +-g input file of Gene region annotation which format like GenePred +-span +-c cluster file input +-o svg output +-h help +USAGE +exit(1); +} \ No newline at end of file