line |
stmt |
bran |
cond |
sub |
pod |
time |
code |
1
|
|
|
|
|
|
|
package App::RL::Command::stat2; |
2
|
13
|
|
|
13
|
|
8803
|
use strict; |
|
13
|
|
|
|
|
32
|
|
|
13
|
|
|
|
|
376
|
|
3
|
13
|
|
|
13
|
|
65
|
use warnings; |
|
13
|
|
|
|
|
31
|
|
|
13
|
|
|
|
|
308
|
|
4
|
13
|
|
|
13
|
|
61
|
use autodie; |
|
13
|
|
|
|
|
27
|
|
|
13
|
|
|
|
|
73
|
|
5
|
|
|
|
|
|
|
|
6
|
13
|
|
|
13
|
|
74953
|
use App::RL -command; |
|
13
|
|
|
|
|
32
|
|
|
13
|
|
|
|
|
129
|
|
7
|
13
|
|
|
13
|
|
4421
|
use App::RL::Common; |
|
13
|
|
|
|
|
30
|
|
|
13
|
|
|
|
|
324
|
|
8
|
|
|
|
|
|
|
|
9
|
13
|
|
|
13
|
|
56
|
use constant abstract => 'coverage statistics on another runlist for runlists'; |
|
13
|
|
|
|
|
25
|
|
|
13
|
|
|
|
|
17154
|
|
10
|
|
|
|
|
|
|
|
11
|
|
|
|
|
|
|
sub opt_spec { |
12
|
|
|
|
|
|
|
return ( |
13
|
8
|
|
|
8
|
1
|
70
|
[ "outfile|o=s", "Output filename. [stdout] for screen" ], |
14
|
|
|
|
|
|
|
[ "op=s", "Operations: intersect, union, diff or xor", { default => "intersect" } ], |
15
|
|
|
|
|
|
|
[ "size|s=s", "chr.sizes", { required => 1 } ], |
16
|
|
|
|
|
|
|
[ "base|b=s", "basename of infile2", ], |
17
|
|
|
|
|
|
|
[ "remove|r", "Remove 'chr0' from chromosome names" ], |
18
|
|
|
|
|
|
|
[ "mk", "First YAML file contains multiple sets of runlists" ], |
19
|
|
|
|
|
|
|
[ "all", "Only write whole genome stats" ], |
20
|
|
|
|
|
|
|
{ show_defaults => 1, } |
21
|
|
|
|
|
|
|
); |
22
|
|
|
|
|
|
|
} |
23
|
|
|
|
|
|
|
|
24
|
|
|
|
|
|
|
sub usage_desc { |
25
|
8
|
|
|
8
|
1
|
52947
|
return "runlist stat2 [options] "; |
26
|
|
|
|
|
|
|
} |
27
|
|
|
|
|
|
|
|
28
|
|
|
|
|
|
|
sub description { |
29
|
1
|
|
|
1
|
1
|
1572
|
my $desc; |
30
|
1
|
|
|
|
|
3
|
$desc .= ucfirst(abstract) . ".\n"; |
31
|
1
|
|
|
|
|
3
|
return $desc; |
32
|
|
|
|
|
|
|
} |
33
|
|
|
|
|
|
|
|
34
|
|
|
|
|
|
|
sub validate_args { |
35
|
6
|
|
|
6
|
1
|
10646
|
my ( $self, $opt, $args ) = @_; |
36
|
|
|
|
|
|
|
|
37
|
6
|
100
|
|
|
|
11
|
if ( @{$args} != 2 ) { |
|
6
|
|
|
|
|
25
|
|
38
|
2
|
|
|
|
|
4
|
my $message = "This command need two input files.\n\tIt found"; |
39
|
2
|
|
|
|
|
3
|
$message .= sprintf " [%s]", $_ for @{$args}; |
|
2
|
|
|
|
|
8
|
|
40
|
2
|
|
|
|
|
5
|
$message .= ".\n"; |
41
|
2
|
|
|
|
|
12
|
$self->usage_error($message); |
42
|
|
|
|
|
|
|
} |
43
|
4
|
|
|
|
|
11
|
for ( @{$args} ) { |
|
4
|
|
|
|
|
12
|
|
44
|
7
|
50
|
|
|
|
265
|
next if lc $_ eq "stdin"; |
45
|
7
|
100
|
|
|
|
29
|
if ( !Path::Tiny::path($_)->is_file ) { |
46
|
1
|
|
|
|
|
114
|
$self->usage_error("The input file [$_] doesn't exist."); |
47
|
|
|
|
|
|
|
} |
48
|
|
|
|
|
|
|
} |
49
|
|
|
|
|
|
|
|
50
|
3
|
50
|
|
|
|
136
|
if ( $opt->{op} =~ /^dif/i ) { |
|
|
50
|
|
|
|
|
|
|
|
50
|
|
|
|
|
|
|
|
0
|
|
|
|
|
|
51
|
0
|
|
|
|
|
0
|
$opt->{op} = 'diff'; |
52
|
|
|
|
|
|
|
} |
53
|
|
|
|
|
|
|
elsif ( $opt->{op} =~ /^uni/i ) { |
54
|
0
|
|
|
|
|
0
|
$opt->{op} = 'union'; |
55
|
|
|
|
|
|
|
} |
56
|
|
|
|
|
|
|
elsif ( $opt->{op} =~ /^int/i ) { |
57
|
3
|
|
|
|
|
9
|
$opt->{op} = 'intersect'; |
58
|
|
|
|
|
|
|
} |
59
|
|
|
|
|
|
|
elsif ( $opt->{op} =~ /^xor/i ) { |
60
|
0
|
|
|
|
|
0
|
$opt->{op} = 'xor'; |
61
|
|
|
|
|
|
|
} |
62
|
|
|
|
|
|
|
else { |
63
|
0
|
|
|
|
|
0
|
Carp::confess "[@{[$opt->{op}]}] invalid\n"; |
|
0
|
|
|
|
|
0
|
|
64
|
|
|
|
|
|
|
} |
65
|
|
|
|
|
|
|
|
66
|
3
|
50
|
|
|
|
13
|
if ( !exists $opt->{base} ) { |
67
|
3
|
|
|
|
|
10
|
$opt->{base} = Path::Tiny::path( $args->[1] )->basename( ".yaml", ".yml" ); |
68
|
|
|
|
|
|
|
} |
69
|
|
|
|
|
|
|
|
70
|
3
|
50
|
|
|
|
288
|
if ( !exists $opt->{outfile} ) { |
71
|
0
|
|
|
|
|
0
|
$opt->{outfile} = Path::Tiny::path( $args->[0] )->absolute . $opt->{op} . ".csv"; |
72
|
|
|
|
|
|
|
} |
73
|
|
|
|
|
|
|
} |
74
|
|
|
|
|
|
|
|
75
|
|
|
|
|
|
|
sub execute { |
76
|
3
|
|
|
3
|
1
|
21
|
my ( $self, $opt, $args ) = @_; |
77
|
|
|
|
|
|
|
|
78
|
|
|
|
|
|
|
#----------------------------# |
79
|
|
|
|
|
|
|
# Loading |
80
|
|
|
|
|
|
|
#----------------------------# |
81
|
3
|
|
|
|
|
24
|
my $length_of = App::RL::Common::read_sizes( $opt->{size}, $opt->{remove} ); |
82
|
|
|
|
|
|
|
|
83
|
|
|
|
|
|
|
# file1 |
84
|
3
|
|
|
|
|
11
|
my $s1_of = {}; |
85
|
3
|
|
|
|
|
7
|
my @keys; |
86
|
3
|
100
|
|
|
|
11
|
if ( $opt->{mk} ) { |
87
|
1
|
|
|
|
|
8
|
my $yml = YAML::Syck::LoadFile( $args->[0] ); |
88
|
1
|
|
|
|
|
146
|
@keys = sort keys %{$yml}; |
|
1
|
|
|
|
|
7
|
|
89
|
|
|
|
|
|
|
|
90
|
1
|
|
|
|
|
4
|
for my $key (@keys) { |
91
|
5
|
|
|
|
|
64
|
$s1_of->{$key} = App::RL::Common::runlist2set( $yml->{$key}, $opt->{remove} ); |
92
|
|
|
|
|
|
|
} |
93
|
|
|
|
|
|
|
} |
94
|
|
|
|
|
|
|
else { |
95
|
2
|
|
|
|
|
5
|
@keys = ("__single"); |
96
|
|
|
|
|
|
|
$s1_of->{__single} |
97
|
2
|
|
|
|
|
13
|
= App::RL::Common::runlist2set( YAML::Syck::LoadFile( $args->[0] ), $opt->{remove} ); |
98
|
|
|
|
|
|
|
} |
99
|
|
|
|
|
|
|
|
100
|
|
|
|
|
|
|
# file2 |
101
|
3
|
|
|
|
|
41
|
my $s2 = App::RL::Common::runlist2set( YAML::Syck::LoadFile( $args->[1] ), $opt->{remove} ); |
102
|
|
|
|
|
|
|
|
103
|
|
|
|
|
|
|
#----------------------------# |
104
|
|
|
|
|
|
|
# Operating |
105
|
|
|
|
|
|
|
#----------------------------# |
106
|
3
|
|
|
|
|
18
|
my $op_result_of = { map { $_ => {} } @keys }; |
|
7
|
|
|
|
|
26
|
|
107
|
|
|
|
|
|
|
|
108
|
3
|
|
|
|
|
10
|
for my $key (@keys) { |
109
|
7
|
|
|
|
|
16
|
my $s1 = $s1_of->{$key}; |
110
|
|
|
|
|
|
|
|
111
|
|
|
|
|
|
|
# give empty set to non-existing chrs |
112
|
7
|
|
|
|
|
16
|
for my $s ( $s1, $s2 ) { |
113
|
14
|
|
|
|
|
63
|
for my $chr ( sort keys %{$length_of} ) { |
|
14
|
|
|
|
|
62
|
|
114
|
114
|
100
|
|
|
|
374
|
if ( !exists $s->{$chr} ) { |
115
|
23
|
|
|
|
|
62
|
$s->{$chr} = App::RL::Common::new_set(); |
116
|
|
|
|
|
|
|
} |
117
|
|
|
|
|
|
|
} |
118
|
|
|
|
|
|
|
} |
119
|
|
|
|
|
|
|
|
120
|
|
|
|
|
|
|
# operate on each chr |
121
|
7
|
|
|
|
|
21
|
for my $chr ( sort keys %{$length_of} ) { |
|
7
|
|
|
|
|
24
|
|
122
|
57
|
|
|
|
|
99
|
my $op = $opt->{op}; |
123
|
57
|
|
|
|
|
157
|
my $op_set = $s1->{$chr}->$op( $s2->{$chr} ); |
124
|
57
|
|
|
|
|
20631
|
$op_result_of->{$key}{$chr} = $op_set; |
125
|
|
|
|
|
|
|
} |
126
|
|
|
|
|
|
|
} |
127
|
|
|
|
|
|
|
|
128
|
|
|
|
|
|
|
# warn YAML::Syck::Dump $s2; |
129
|
|
|
|
|
|
|
# warn YAML::Syck::Dump $op_result_of; |
130
|
|
|
|
|
|
|
|
131
|
|
|
|
|
|
|
#----------------------------# |
132
|
|
|
|
|
|
|
# Calcing |
133
|
|
|
|
|
|
|
#----------------------------# |
134
|
3
|
|
|
|
|
6
|
my $out_fh; |
135
|
3
|
50
|
|
|
|
17
|
if ( lc( $opt->{outfile} ) eq "stdout" ) { |
136
|
3
|
|
|
|
|
15
|
$out_fh = *STDOUT; |
137
|
|
|
|
|
|
|
} |
138
|
|
|
|
|
|
|
else { |
139
|
0
|
|
|
|
|
0
|
open $out_fh, ">", $opt->{outfile}; |
140
|
|
|
|
|
|
|
} |
141
|
|
|
|
|
|
|
|
142
|
|
|
|
|
|
|
my $header = sprintf "key,chr,chr_length,size,%s_length,%s_size,c1,c2,ratio\n", |
143
|
3
|
|
|
|
|
26
|
$opt->{base}, $opt->{base}; |
144
|
3
|
100
|
|
|
|
10
|
if ( $opt->{mk} ) { |
145
|
1
|
50
|
|
|
|
5
|
if ( $opt->{all} ) { |
146
|
1
|
|
|
|
|
5
|
$header =~ s/chr\,//; |
147
|
|
|
|
|
|
|
} |
148
|
1
|
|
|
|
|
3
|
my @lines = ($header); |
149
|
|
|
|
|
|
|
|
150
|
1
|
|
|
|
|
3
|
for my $key (@keys) { |
151
|
|
|
|
|
|
|
my @key_lines |
152
|
5
|
|
|
|
|
16
|
= csv_lines( $s1_of->{$key}, $length_of, $s2, $op_result_of->{$key}, $opt->{all} ); |
153
|
5
|
|
|
|
|
15
|
$_ = "$key,$_" for @key_lines; |
154
|
5
|
|
|
|
|
11
|
push @lines, @key_lines; |
155
|
|
|
|
|
|
|
} |
156
|
|
|
|
|
|
|
|
157
|
1
|
|
|
|
|
4
|
print {$out_fh} $_ for @lines; |
|
6
|
|
|
|
|
78
|
|
158
|
|
|
|
|
|
|
} |
159
|
|
|
|
|
|
|
else { |
160
|
2
|
|
|
|
|
14
|
$header =~ s/key\,//; |
161
|
2
|
100
|
|
|
|
12
|
if ( $opt->{all} ) { |
162
|
1
|
|
|
|
|
5
|
$header =~ s/chr\,//; |
163
|
|
|
|
|
|
|
} |
164
|
2
|
|
|
|
|
6
|
my @lines = ($header); |
165
|
|
|
|
|
|
|
|
166
|
|
|
|
|
|
|
push @lines, |
167
|
|
|
|
|
|
|
csv_lines( $s1_of->{__single}, $length_of, $s2, $op_result_of->{__single}, |
168
|
2
|
|
|
|
|
11
|
$opt->{all} ); |
169
|
|
|
|
|
|
|
|
170
|
2
|
|
|
|
|
7
|
print {$out_fh} $_ for @lines; |
|
20
|
|
|
|
|
263
|
|
171
|
|
|
|
|
|
|
} |
172
|
|
|
|
|
|
|
|
173
|
3
|
|
|
|
|
48
|
close $out_fh; |
174
|
|
|
|
|
|
|
} |
175
|
|
|
|
|
|
|
|
176
|
|
|
|
|
|
|
sub csv_lines { |
177
|
7
|
|
|
7
|
0
|
10
|
my $set_of = shift; |
178
|
7
|
|
|
|
|
10
|
my $length_of = shift; |
179
|
7
|
|
|
|
|
15
|
my $s2 = shift; |
180
|
7
|
|
|
|
|
9
|
my $op_result_of = shift; |
181
|
7
|
|
|
|
|
12
|
my $all = shift; |
182
|
|
|
|
|
|
|
|
183
|
7
|
|
|
|
|
12
|
my @lines; |
184
|
7
|
|
|
|
|
16
|
my ( $all_length, $all_size, $all_s2_length, $all_s2_size, ); |
185
|
7
|
|
|
|
|
13
|
for my $chr ( sort keys %{$length_of} ) { |
|
7
|
|
|
|
|
33
|
|
186
|
57
|
|
|
|
|
95
|
my $length = $length_of->{$chr}; |
187
|
57
|
|
|
|
|
138
|
my $size = $set_of->{$chr}->size; |
188
|
|
|
|
|
|
|
|
189
|
57
|
|
|
|
|
10076
|
my $s2_length = $s2->{$chr}->size; |
190
|
57
|
|
|
|
|
1275
|
my $s2_size = $op_result_of->{$chr}->size; |
191
|
|
|
|
|
|
|
|
192
|
57
|
|
|
|
|
792
|
$all_length += $length; |
193
|
57
|
|
|
|
|
78
|
$all_size += $size; |
194
|
57
|
|
|
|
|
70
|
$all_s2_length += $s2_length; |
195
|
57
|
|
|
|
|
78
|
$all_s2_size += $s2_size; |
196
|
|
|
|
|
|
|
|
197
|
57
|
|
|
|
|
80
|
my $c1 = $size / $length; |
198
|
57
|
100
|
|
|
|
110
|
my $c2 = $s2_length == 0 ? 0 : $s2_size / $s2_length; |
199
|
57
|
100
|
|
|
|
108
|
my $ratio = $c1 == 0 ? 0 : $c2 / $c1; |
200
|
|
|
|
|
|
|
|
201
|
57
|
|
|
|
|
399
|
my $line = sprintf "%s,%d,%d,%d,%d,%.4f,%.4f,%.4f\n", $chr, $length, $size, $s2_length, |
202
|
|
|
|
|
|
|
$s2_size, $c1, $c2, $ratio; |
203
|
57
|
|
|
|
|
125
|
push @lines, $line; |
204
|
|
|
|
|
|
|
} |
205
|
|
|
|
|
|
|
|
206
|
7
|
|
|
|
|
15
|
my $all_c1 = $all_size / $all_length; |
207
|
7
|
50
|
|
|
|
24
|
my $all_c2 = $all_s2_length == 0 ? 0 : $all_s2_size / $all_s2_length; |
208
|
7
|
50
|
|
|
|
18
|
my $all_ratio = $all_c1 == 0 ? 0 : $all_c2 / $all_c1; |
209
|
|
|
|
|
|
|
|
210
|
|
|
|
|
|
|
# only keep whole genome |
211
|
7
|
|
|
|
|
40
|
my $all_line = sprintf "all,%d,%d,%d,%d,%.4f,%.4f,%.4f\n", $all_length, $all_size, |
212
|
|
|
|
|
|
|
$all_s2_length, $all_s2_size, $all_c1, $all_c2, $all_ratio; |
213
|
7
|
100
|
|
|
|
20
|
if ($all) { |
214
|
6
|
|
|
|
|
13
|
@lines = (); |
215
|
6
|
|
|
|
|
22
|
$all_line =~ s/all,//; |
216
|
|
|
|
|
|
|
} |
217
|
7
|
|
|
|
|
15
|
push @lines, $all_line; |
218
|
|
|
|
|
|
|
|
219
|
7
|
|
|
|
|
24
|
return @lines; |
220
|
|
|
|
|
|
|
} |
221
|
|
|
|
|
|
|
|
222
|
|
|
|
|
|
|
1; |