-
Notifications
You must be signed in to change notification settings - Fork 1
/
count_every_base_in_bam_compare_output.pl~
executable file
·70 lines (67 loc) · 1.76 KB
/
count_every_base_in_bam_compare_output.pl~
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
#!/usr/bin/perl -w
use strict;
use Data::Dumper;
my @args = @ARGV;
if (!defined @ARGV ){
die "script.pl file id pingCount file id pingCount file id pingCount\n";
}
my %args;
for (my $i =0; $i < @args ; $i+=3){
$args{$args[$i]}{id}=$args[$i+1];
$args{$args[$i]}{ping}=$args[$i+2];
}
my %pos;
foreach my $file (keys %args){
open IN, $file;
my $id = $args{$file}{id};
my $ping = $args{$file}{ping};
while (my $line = <IN>){
chomp $line;
my $diff=0;
my ($pos, $A, $T, $G, $C, $N, $total)= split /\t/, $line;
foreach my $nt_count ($A, $T, $G, $C ){
my ($nt , $count) = split /\:/ , $nt_count;
#if ($count > ($total/$ping * .2)){
if ($count > ($total/$ping * 0.5)){
$diff++;
}
if ($nt eq 'A'){
$A = $count;
}if ($nt eq 'T'){
$T = $count;
}if ($nt eq 'G'){
$G = $count;
}if ($nt eq 'C'){
$C = $count;
}
}
if ($diff > 1){
$pos{$pos}{$id}{A}=$A;
$pos{$pos}{$id}{T}=$T;
$pos{$pos}{$id}{G}=$G;
$pos{$pos}{$id}{C}=$C;
$pos{$pos}{$id}{total}=$total;
}
}
}
foreach my $pos (sort { (split /\(/, $a)[0] <=> (split /\(/ ,$b)[0]} keys %pos){
my $line = "$pos ::";
my $diffcount = 0;
my $HEG4_count= 0;
foreach my $id (sort keys %{$pos{$pos}}){
$HEG4_count++ if $id =~ /HEG4/;
my $total = $pos{$pos}{$id}{total};
$line .= "\n\t$id($total)\t";
$diffcount++;
foreach my $nt (sort keys %{$pos{$pos}{$id}}){
next if $nt eq 'total';
my $count = $pos{$pos}{$id}{$nt};
$line .= "$nt:$count\t";
}
}
$line =~s/\t$/\n/;
$line =~s/\s::/:($diffcount)/;
next if $HEG4_count == 1;
next if $line =~ /HEG4_0/ and $diffcount eq 1;
print $line;
}