The Perl and Raku Conference 2025: Greenville, South Carolina - June 27-29 Learn more

MNSMAR / CLIPSeqTools-1.0.0 / lib / CLIPSeqTools / PlotApp / genomic_distribution.pm 

            

              
—
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125

              
=head1 NAME
CLIPSeqTools::PlotApp::genomic_distribution - Create plots for script genomic_distribution.
=head1 SYNOPSIS
clipseqtools-plot genomic_distribution [options/parameters]
=head1 DESCRIPTION
Create plots for script genomic_distribution.
=head1 OPTIONS
  Input.
    --file <Str>           input file.
  Output
    --o_prefix <Str>       output path prefix. Script will create and add
                           extension to path. [Default: ./]
    -v --verbose           print progress lines and extra information.
    -h -? --usage --help   print help message
=cut
package CLIPSeqTools::PlotApp::genomic_distribution;
$CLIPSeqTools::PlotApp::genomic_distribution::VERSION = '1.0.0';
# Make it an app command
use MooseX::App::Command;
extends 'CLIPSeqTools::PlotApp';
#######################################################################
#######################   Load External modules   #####################
#######################################################################
use Modern::Perl;
use autodie;
use namespace::autoclean;
use Statistics::R;
#######################################################################
#######################   Command line options   ######################
#######################################################################
option 'file' => (
        is            => 'rw',
        isa           => 'Str',
        required      => 1,
        documentation => 'input file.',
);
#######################################################################
##########################   Consume Roles   ##########################
#######################################################################
with
        "CLIPSeqTools::Role::Option::OutputPrefix" => {
                -alias    => { validate_args => '_validate_args_for_output_prefix' },
                -excludes => 'validate_args',
        };
#######################################################################
########################   Interface Methods   ########################
#######################################################################
sub validate_args {
        my ($self) = @_;
        $self->_validate_args_for_output_prefix;
}
sub run {
        my ($self) = @_;
        warn "Validating arguments\n" if $self->verbose;
        $self->validate_args();
        warn "Creating output path\n" if $self->verbose;
        $self->make_path_for_output_prefix();
        warn "Creating plots with R\n" if $self->verbose;
        $self->run_R;
}
sub run_R {
        my ($self) = @_;
        my $figfile = $self->o_prefix . 'genomic_distribution.pdf';
        # Start R
        my $R = Statistics::R->new();
        # Pass arguments to R
        $R->set('ifile', $self->file);
        $R->set('figfile', $figfile);
        # Read table with data - Do exra calulations
        $R->run(q{idata = read.delim(ifile)});
        $R->run(q{idata$percent = (idata$count / idata$total) * 100});
        # Do plots
        $R->run(q{pdf(figfile, width=14)});
        $R->run(q{par(mfrow = c(1, 2), mar=c(10, 4.1, 3.1, 2.1));});
        $R->run(q{mypal = c(rep("black",2), rep("darkgrey",3), rep("grey",3),
                rep("lightgrey",3), rep("lightblue",3))});
        $R->run(q{x = barplot(idata$percent, names.arg=idata$category, col=mypal,
                ylim=c(0,100), ylab="Percent of total reads", xaxt="n");});
        $R->run(q{text(x=x, y=-max(idata$percent)/100, idata$category, xpd=TRUE,
                srt=50, adj=c(1, 1))});
        $R->run(q{x = barplot(height=idata$count, names.arg=idata$category,
                col=mypal, ylab="Number of reads", xaxt="n");});
        $R->run(q{text(x=x, y=-max(idata$count)/100, idata$category, xpd=TRUE,
                srt=50, adj=c(1, 1))});
        $R->run(q{graphics.off()});
        # Close R
        $R->stop();
}
1;