The Perl Toolchain Summit 2025 Needs You: You can help 🙏 Learn more

MNSMAR / CLIPSeqTools-0.1.0-TRIAL / lib / CLIPSeqTools / App / size_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
126
127
128
129
130
131
132
133
134
135

              
=head1 NAME
CLIPSeqTools::App::size_distribution - Measure size distribution for reads.
=head1 SYNOPSIS
clipseqtools size_distribution [options/parameters]
=head1 DESCRIPTION
Measure size distribution of aligned reads.
=head1 OPTIONS
  Input options for library.
    --driver <Str>         driver for database connection (eg. mysql,
                           SQLite).
    --database <Str>       database name or path to database file for file
                           based databases (eg. SQLite).
    --table <Str>          database table.
    --host <Str>           hostname for database connection.
    --user <Str>           username for database connection.
    --password <Str>       password for database connection.
    --records_class <Str>  type of records stored in database.
    --filter <Filter>      filter library. May be used multiple times.
                           Syntax: column_name="pattern"
                           e.g. keep reads with deletions AND not repeat
                                masked AND longer than 31
                                --filter deletion="def" 
                                --filter rmsk="undef" .
                                --filter query_length=">31".
                           Operators: >, >=, <, <=, =, !=, def, undef
  Output
    --o_prefix <Str>       output path prefix. Script will create and add
                           extension to path. [Default: ./]
  Other options.
    --plot                 call plotting script to create plots.
    -v --verbose           print progress lines and extra information.
    -h -? --usage --help   print help message
=cut
package CLIPSeqTools::App::size_distribution;
$CLIPSeqTools::App::size_distribution::VERSION = '0.1.0'; # TRIAL
# Make it an app command
use MooseX::App::Command;
extends 'CLIPSeqTools::App';
#######################################################################
#######################   Load External modules   #####################
#######################################################################
use Modern::Perl;
use autodie;
use namespace::autoclean;
use List::Util qw(min max);
#######################################################################
##########################   Consume Roles   ##########################
#######################################################################
with 
        "CLIPSeqTools::Role::Option::Library" => {
                -alias    => { validate_args => '_validate_args_for_library' },
                -excludes => 'validate_args',
        },
        "CLIPSeqTools::Role::Option::Plot" => {
                -alias    => { validate_args => '_validate_args_for_plot' },
                -excludes => 'validate_args',
        },
        "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_library;
        $self->_validate_args_for_plot;
        $self->_validate_args_for_output_prefix;
}
sub run {
        my ($self) = @_;
         
        warn "Starting analysis: size_distribution\n";
         
        warn "Validating arguments\n" if $self->verbose;
        $self->validate_args();
         
        warn "Creating reads collection\n" if $self->verbose;
        my $reads_collection = $self->reads_collection;
        $reads_collection->schema->storage->debug(1) if $self->verbose > 1;
         
        warn "Measuring aligned read sizes\n" if $self->verbose;
        my %size_count;
        $reads_collection->foreach_record_do( sub {
                my ($rec) = @_;
                 
                return 0 if $rec->cigar =~ /N/; # throw away reads with huge gaps (introns)
                 
                $size_count{$rec->length} += $rec->copy_number;
                 
                return 0;
        });
        my $min_length = min(keys %size_count);
        my $max_length = max(keys %size_count);
        warn "Creating output path\n" if $self->verbose;
        $self->make_path_for_output_prefix();
         
        warn "Printing results\n" if $self->verbose;
        open (my $OUT1, '>', $self->o_prefix.'size_distribution.tab');
        say $OUT1 join("\t", 'size', 'count');
        say $OUT1 join("\t", $_, $size_count{$_} || 0) for ($min_length..$max_length);
        close $OUT1;
         
        if ($self->plot) {
                warn "Creating plot\n" if $self->verbose;
                CLIPSeqTools::PlotApp->initialize_command_class('CLIPSeqTools::PlotApp::size_distribution', 
                        file     => $self->o_prefix.'size_distribution.tab',
                        o_prefix => $self->o_prefix
                )->run();
        }
}
1;