Mercurial > repos > romaingred > pirna_pipeline
comparison bin/Rcall.pm @ 0:198009598544 draft
Uploaded
author | romaingred |
---|---|
date | Wed, 11 Oct 2017 09:57:58 -0400 |
parents | |
children | ee6b4b2072a9 |
comparison
equal
deleted
inserted
replaced
-1:000000000000 | 0:198009598544 |
---|---|
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 = $size_hashR->{$k} * 100 / $size; | |
19 $abs[$i] = $k ; $ord[$i] = $percentage; $i++; | |
20 } | |
21 my $abs = join (",", @abs ); | |
22 my $ord = join (",", @ord ); | |
23 if (scalar(@abs) != 0) | |
24 { | |
25 | |
26 my $R = Statistics::R->new(); | |
27 $R->startR; | |
28 $R->send( | |
29 qq`library(ggplot2) | |
30 percentage = c($ord) | |
31 size =c($abs) | |
32 min = min(size) | |
33 max = max(size) | |
34 dat = data.frame(size,percentage) | |
35 png(filename=\"$out_png\", width = 800, height = 480) | |
36 c = ggplot(dat,aes(size,percentage)) | |
37 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") ) | |
38 dev.off()`); | |
39 $R->stopR(); | |
40 | |
41 } | |
42 } | |
43 | |
44 sub bg_to_png | |
45 { | |
46 my ( $fai, $bgP, $bgM, $dir, $sb ) = @_; | |
47 my $R = Statistics::R->new(); | |
48 $R->startR; | |
49 $R->send( | |
50 qq`library('Sushi') | |
51 fai =read.table("$fai") | |
52 bgP = read.table("$bgP") | |
53 bgM = read.table("$bgM") | |
54 fun = function(chr,end) { | |
55 jpeg( paste0("$dir",as.character(chr),".png"), quality=100) | |
56 par(mfrow=c(2,1),mar=c(1,10,1,3)) | |
57 plotBedgraph(bgP, chrom=chr,chromstart=0,chromend=end,transparency=.50, color=SushiColors(2)(2)[1]) | |
58 axis(side=2,las=2,tcl=.2) | |
59 mtext("Scaled Read Depth",side=2,line=4,cex=1,font=2) | |
60 plotBedgraph(bgM, chrom=chr,chromstart=0,chromend=end,transparency=.50, flip=TRUE, color=SushiColors(2)(2)[2]) | |
61 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) | |
62 axis(side=2,las=2,tcl=.2,at=pretty(par("yaxp")[c(1,2)]),labels=-1*pretty(par("yaxp")[c(1,2)])) | |
63 mtext("Scaled Read Depth",side=2,line=4.5,cex=1,font=2) | |
64 dev.off() | |
65 } | |
66 mapply( fun, fai\$V1, fai\$V2)`); | |
67 $R->stopR(); | |
68 } | |
69 | |
70 sub pie_chart | |
71 { | |
72 my $dir = shift; | |
73 my $in = $dir.'repartition.txt'; | |
74 my $out = $dir.'pie_chart.png'; | |
75 | |
76 my $R = Statistics::R->new(); | |
77 $R->startR; | |
78 $R->send( | |
79 qq` | |
80 library(plotrix) | |
81 library(RColorBrewer) | |
82 R =read.table("$in",header=T) | |
83 values = round(R\$percentage) | |
84 keys = R\$type | |
85 lab = paste(values, "%", sep="") | |
86 png("$out") | |
87 colors <- brewer.pal(7,"Paired") | |
88 pie(values, col=colors, labels=lab, clockwise=TRUE) | |
89 legend("bottom", legend = keys, fill=colors, bty="n", ncol = 3) | |
90 par(mai = c(0,0,0,0)) | |
91 layout(c(1,2),heights=c(0.3,1)) | |
92 plot.new() | |
93 legend("bottom", legend = keys, fill=colors, bty="n",ncol = 3) | |
94 pie(values, col=colors, labels=lab, clockwise=TRUE) | |
95 dev.off()` | |
96 ); | |
97 $R->stopR(); | |
98 } | |
99 | |
100 1; |