Mercurial > repos > romaingred > srnapipe
comparison bin/Rcall.pm @ 4:bf3c8bf8b819 draft
Deleted selected files
| author | romaingred |
|---|---|
| date | Wed, 13 Dec 2017 10:24:17 -0500 |
| parents | |
| children |
comparison
equal
deleted
inserted
replaced
| 3:e9c9cee888cf | 4:bf3c8bf8b819 |
|---|---|
| 1 package Rcall; | |
| 2 | |
| 3 use strict; | |
| 4 use warnings; | |
| 5 use Statistics::R; | |
| 6 | |
| 7 use Exporter; | |
| 8 our @ISA = qw(Exporter); | |
| 9 our @EXPORT_OK = qw( &histogram &pie_chart &bg_to_png ); | |
| 10 | |
| 11 sub histogram | |
| 12 { | |
| 13 my ($size_hashR, $out_png, $size) = @_; | |
| 14 my (@abs, @ord); | |
| 15 my $i = 0; | |
| 16 foreach my $k (sort {$a <=> $b} keys %{$size_hashR}) | |
| 17 { | |
| 18 my $percentage = 0; | |
| 19 $percentage = $size_hashR->{$k} * 100 / $size if $size != 0; | |
| 20 $abs[$i] = $k ; $ord[$i] = $percentage; $i++; | |
| 21 } | |
| 22 my $abs = join (",", @abs ); | |
| 23 my $ord = join (",", @ord ); | |
| 24 if (scalar(@abs) != 0) | |
| 25 { | |
| 26 | |
| 27 my $R = Statistics::R->new(); | |
| 28 $R->startR; | |
| 29 $R->send( | |
| 30 qq`library(ggplot2) | |
| 31 percentage = c($ord) | |
| 32 size =c($abs) | |
| 33 min = min(size) | |
| 34 max = max(size) | |
| 35 dat = data.frame(size,percentage) | |
| 36 png(filename=\"$out_png\", width = 800, height = 480) | |
| 37 c = ggplot(dat,aes(size,percentage)) | |
| 38 c + geom_bar(stat="identity") + scale_x_continuous(breaks=min:max)+theme( axis.text.x = element_text(angle=90, hjust=0.5, size=20), axis.text.y = element_text( size=20 ), axis.title.x = element_text( size=25, face="bold"), axis.title.y = element_text( size=25, face="bold") ) | |
| 39 dev.off()`); | |
| 40 $R->stopR(); | |
| 41 | |
| 42 } | |
| 43 } | |
| 44 | |
| 45 sub bg_to_png | |
| 46 { | |
| 47 my ( $fai, $bgP, $bgM, $dir, $sb ) = @_; | |
| 48 my $R = Statistics::R->new(); | |
| 49 $R->startR; | |
| 50 $R->send( | |
| 51 qq`library('Sushi') | |
| 52 fai =read.table("$fai") | |
| 53 if ( file.info("$bgP")\$size !=0 ) | |
| 54 { | |
| 55 bgP = read.table("$bgP") | |
| 56 } else { bgP = data.frame(factor(),integer()) } | |
| 57 | |
| 58 if ( file.info("$bgM")\$size !=0 ) | |
| 59 { | |
| 60 bgM = read.table("$bgM") | |
| 61 } else { bgM = data.frame(factor(),integer()) } | |
| 62 | |
| 63 f_both = function(chr,end) { | |
| 64 jpeg( paste0("$dir",as.character(chr),".png"), quality=100) | |
| 65 par(mfrow=c(2,1),mar=c(1,10,1,3)) | |
| 66 plotBedgraph(bgP, chrom=chr,chromstart=0,chromend=end,transparency=.50, color=SushiColors(2)(2)[1]) | |
| 67 axis(side=2,las=2,tcl=.2) | |
| 68 mtext("Scaled Read Depth",side=2,line=4,cex=1,font=2) | |
| 69 plotBedgraph(bgM, chrom=chr,chromstart=0,chromend=end,transparency=.50, flip=TRUE, color=SushiColors(2)(2)[2]) | |
| 70 labelgenome(chrom=chr,chromstart=0,chromend=end,side=3,n=3,scale="$sb", line=0, chromline = 0.5, scaleline = 0.5, scaleadjust =1.05, chromadjust = -0.4) | |
| 71 axis(side=2,las=2,tcl=.2,at=pretty(par("yaxp")[c(1,2)]),labels=-1*pretty(par("yaxp")[c(1,2)])) | |
| 72 mtext("Scaled Read Depth",side=2,line=4.5,cex=1,font=2) | |
| 73 dev.off() | |
| 74 } | |
| 75 | |
| 76 f_plus = function(chr,end) { | |
| 77 jpeg( paste0("$dir",as.character(chr),".png"), quality=100) | |
| 78 plotBedgraph(bgP, chrom=chr,chromstart=0,chromend=end,transparency=.50, color=SushiColors(2)(2)[1]) | |
| 79 labelgenome(chrom=chr,chromstart=0,chromend=end,n=3,scale="$sb", line=0, chromline = 0.5, scaleline = 0.5, scaleadjust =1.05, chromadjust = -0.4) | |
| 80 axis(side=2,las=2,tcl=.2) | |
| 81 mtext("Scaled Read Depth",side=2,line=4,cex=1,font=2) | |
| 82 dev.off() | |
| 83 } | |
| 84 | |
| 85 f_minus = function(chr,end) { | |
| 86 jpeg( paste0("$dir",as.character(chr),".png"), quality=100) | |
| 87 plotBedgraph(bgM, chrom=chr,chromstart=0,chromend=end,transparency=.50, flip=TRUE, color=SushiColors(2)(2)[2]) | |
| 88 labelgenome(chrom=chr,chromstart=0,chromend=end,n=3,scale="$sb", line=0, chromline = 0.5, scaleline = 0.5, scaleadjust =1.05, chromadjust = -0.4) | |
| 89 axis(side=2,las=2,tcl=.2,at=pretty(par("yaxp")[c(1,2)]),labels=-1*pretty(par("yaxp")[c(1,2)])) | |
| 90 mtext("Scaled Read Depth",side=2,line=4.5,cex=1,font=2) | |
| 91 dev.off() | |
| 92 } | |
| 93 | |
| 94 fai_b = fai[fai\$V1 %in% intersect(bgM\$V1,bgP\$V1), ] | |
| 95 mapply( f_both, fai_b\$V1, fai_b\$V2) | |
| 96 | |
| 97 fai_p = fai[fai\$V1 %in% setdiff(bgP\$V1,bgM\$V1), ] | |
| 98 mapply( f_plus, fai_p\$V1, fai_p\$V2) | |
| 99 | |
| 100 fai_m = fai[fai\$V1 %in% setdiff(bgM\$V1,bgP\$V1), ] | |
| 101 mapply( f_minus, fai_m\$V1, fai_m\$V2) `); | |
| 102 | |
| 103 $R->stopR(); | |
| 104 } | |
| 105 | |
| 106 sub pie_chart | |
| 107 { | |
| 108 my $dir = shift; | |
| 109 my $in = $dir.'repartition.txt'; | |
| 110 my $out = $dir.'pie_chart.png'; | |
| 111 | |
| 112 my $R = Statistics::R->new(); | |
| 113 $R->startR; | |
| 114 $R->send( | |
| 115 qq` | |
| 116 library(plotrix) | |
| 117 library(RColorBrewer) | |
| 118 R =read.table("$in",header=T) | |
| 119 values = round(R\$percentage) | |
| 120 keys = R\$type | |
| 121 lab = paste(values, "%", sep="") | |
| 122 png("$out") | |
| 123 colors <- brewer.pal(7,"Paired") | |
| 124 pie(values, col=colors, labels=lab, clockwise=TRUE) | |
| 125 legend("bottom", legend = keys, fill=colors, bty="n", ncol = 3) | |
| 126 par(mai = c(0,0,0,0)) | |
| 127 layout(c(1,2),heights=c(0.3,1)) | |
| 128 plot.new() | |
| 129 legend("bottom", legend = keys, fill=colors, bty="n",ncol = 3) | |
| 130 pie(values, col=colors, labels=lab, clockwise=TRUE) | |
| 131 dev.off()` | |
| 132 ); | |
| 133 $R->stopR(); | |
| 134 } | |
| 135 | |
| 136 1; |
