求指导

那只猴子

有如下格式文件，类型bed文件：
AF043090.1      0       2790    2790
AF192261.1      5       203     1054
AF192261.1      301     1028    1054
AF253472.1      0       4175    4309
AF439974.1      16      106     1409
AF439974.1      248     641     1409
AF439974.1      694     1408    1409
AF460219.1      0       2374    4098
AF460219.1      2383    2594    4098
AF460219.1      2669    2759    4098

第一列是ID；
第二列是起始位点，从0开始；
第三列是终止位点；
第四是该ID的全长；


我想利用perl程序，得到不在该文件记录的位置信息，即想得到以下结果：
AF192261.1 0-4,204-300，1029-1054
AF253472.1 4176-4309
AF439974.1 0-15,107-247,642-693
AF460219.1 2375-2382，2595-2678，2760-4098

请大家指教。

huang6894

回复 1# 那只猴子


    菜鸟写法，仅作参考：

1. #!/usr/bin/perl -w
   
2. use 5.010;
   
3. my (%le,%t);
   
4. while(<DATA>){
   
5.                 chomp;
   
6.                 my($id, $s, $e, $l) = (split);
   
7.                 $le{$id} = $l;
   
8.                 push @{$t{$id}},[$s, $e];
   
9. }
   
10. foreach my $id(keys %le){
    
11.                 my @tmp = sort { $a->[0] <=> $b->[0] }@{$t{$id}};
    
12.                 my @all;
    
13.                 for my $r (@tmp) {
    
14.                         last if $r->[1] - $r->[0] > $le{$id};
    
15.                         for my $i(0..$le{$id}){
    
16.                         push @all,$i unless (($r->[0] <= $i) and  ($i <= $r->[1]));
    
17.                         }
    
18.                 }
    
19.                 my $s = join(",",@all);
    
20.                 $s ||= "NA";
    
21.                 1 while $s =~ s/(-)?(?<!\d)(\d+),(\d++)(?(?{$3 != $2+1})(*F))/defined $1 ?  "-$3" : "$2-$3" /e;
    
22.                 say "$id\t$s";
    
23. }
    
24. 
    
25. __DATA__
    
26. AF043090.1      0       2790    2790
    
27. AF192261.1      5       203     1054
    
28. AF192261.1      301     1028    1054
    
29. AF253472.1      0       4175    4309
    
30. AF439974.1      16      106     1409
    
31. AF439974.1      248     641     1409
    
32. AF439974.1      694     1408    1409
    
33. AF460219.1      0       2374    4098
    
34. AF460219.1      2383    2594    4098
    
35. AF460219.1      2669    2759    4098
    

复制代码

----------------------------额，有bug，不改了，看三楼吧------------------------

qw0148

看看行不行，有个bug是会打印1409-1409

#!/usr/bin/perl
use strict;
my %hash;
my %hash1;
while(<DATA>){
                chomp;
                my ($id,$beg,$end,$last)=split /\s+/;
                $hash{$id}{$beg}=$end;
                $hash1{$id}=$last;
}
foreach my $ID(keys %hash){
                my $first=-1;
                next if($hash{$ID}{0}==$hash1{$ID}); #全部记录则不打印ID;
                print "$ID " ;
                for my $BEG(sort {$a<=>$b} keys %{$hash{$ID}}){
                                        printf "%d-%d,",$first+1,$BEG-1 unless($BEG==0);
                                        $first=$hash{$ID}{$BEG};
                                }
                        ($first < $hash1{$ID})?printf "%d-%d\n",$first+1,$hash1{$ID}:print "\n";                       
}
__DATA__
AF043090.1      0       2790    2790
AF192261.1      5       203     1054
AF192261.1      301     1028    1054
AF253472.1      0       4175    4309
AF439974.1      16      106     1409
AF439974.1      248     641     1409
AF439974.1      694     1408    1409
AF460219.1      0       2374    4098
AF460219.1      2383    2594    4098
AF460219.1      2669    2759    4098


结果是
AF253472.1 4176-4309
AF439974.1 0-15,107-247,642-693,1409-1409
AF192261.1 0-4,204-300,1029-1054
AF460219.1 2375-2382,2595-2668,2760-4098