- 论坛徽章:
- 32
|
本帖最后由 yestreenstars 于 2014-08-19 12:28 编辑
Debug by yourself.- #!/usr/bin/perl
- use warnings;
- use strict;
- my %hash1;
- open(FILE1, "b") or die;
- push(@{$hash1{(split)[1]}{(split)[4]}{(split)[5]}}, [split]) while <FILE1>;
- close(FILE1);
- my($chr, %hash2, @result);
- open(FILE2, "a") or die;
- while(<FILE2>){
- my @field = split;
- $chr = $field[1] and $hash2{$chr} = join("\t", (split)[0, 8]) and next if $#field == 10;
- next if /^\s*$/;
- for my $i (keys %{$hash1{$chr}}){
- for my $j (keys %{$hash1{$chr}{$i}}){
- my @num;
- while(@{$hash1{$chr}{$i}{$j}}){
- my $temp = $hash1{$chr}{$i}{$j}[0][2];
- shift @{$hash1{$chr}{$i}{$j}} and next if $temp < $field[1];
- last if $temp > $field[2];
- push(@num, $temp);
- shift @{$hash1{$chr}{$i}{$j}};
- }
- push(@result, [16642, $chr, @num[0, -1], $i, $j, $hash2{$chr}, $field[0]]) if @num && ($num[-1] - $num[0] + 1) / ($field[2] - $field[1] + 1) > 0.66;
- }
- }
- }
- close(FILE2);
- while(@result){
- if($result[1]){
- if(join("\t", @{$result[0]}[1, 4, 5]) eq join("\t", @{$result[1]}[1, 4, 5])){
- if($result[0][3] + 1 == $result[1][2]){
- $result[1][2] = $result[0][2];
- $result[1][-1] = "$result[0][-1]:$result[1][-1]";
- }else{
- print join("\t", join("\t", @{$result[0]}[0, 1]), "$result[0][2]-$result[0][3]", join("\t", @{$result[0]}[4..7])), $/;
- }
- shift @result;
- next;
- }
- }
- print join("\t", join("\t", @{$result[0]}[0, 1]), "$result[0][2]-$result[0][3]", join("\t", @{$result[0]}[4..7])), $/;
- shift @result;
- }
复制代码- [root@localhost ~]# cat a
- NM_001145717 chr6 + 36238237 36276372 36238237 36275493 8 PNPLA1 cmpl cmpl
- CDS1 36238437 36238441 0 1 1
- intron1 36238442 36238446 . 205 206
- CDS2 36259447 36259329 1 2 206
- intron2 36259330 36260837 . 438 439
- CDS3 36260838 36260903 0 3 439
- intron3 36260904 36261966 . 504 505
- CDS4 36261967 36262176 0 4 505
- intron4 36262177 36263140 . 714 715
- CDS5 36263141 36263201 0 5 715
- intron5 36263202 36269637 . 775 776
- CDS6 36269638 36270246 1 6 776
- intron6 36270247 36274068 . 1384 1385
- CDS7 36274069 36274153 1 7 1385
- intron7 36274154 36275363 . 1469 1470
- CDS8 36275364 36275493 2 8 1470
- 3-UTR 36275494 36276372 . 8 *1 *879
- NM_002364 chrX + 2336324 2336356 2336328 2336357 2 MAGEB2 cmpl cmpl
- 5-UTR 2336324 2336326 . 1 -103 -6
- intron1 2336327 2336329 . -6 -5
- 5-UTR 2336330 2336337 . 2 -5 -1
- CDS1 2336338 2336347 0 2 1
- 3-UTR 2336348 2336356 . 2 *1 *549
- NM_019056 chrX - 47001615 47004609 47001716 47004078 3 NDUFB11 cmpl cmpl
- 5-UTR 47004079 47004609 . 1 -1 -531
- CDS1 47003872 47004078 0 1 1
- intron1 47002144 47003871 . 208 207
- CDS2 47001983 47002143 0 2 208
- intron2 47001840 47001982 . 369 368
- CDS3 47001716 47001839 2 3 369
- 3-UTR 47001615 47001715 . 3 *101 *1
- [root@localhost ~]# cat b
- 16642 chr6 36238435 18 HOM del
- 16642 chr6 36238437 17 HOM del
- 16642 chr6 36238438 18 HOM del
- 16642 chr6 36238439 17 HOM del
- 16642 chr6 36238440 17 HOM del
- 16642 chr6 36238441 24 HOM del
- 16642 chr6 36238442 24 HOM del
- 16642 chr6 36238443 24 HOM del
- 16642 chr6 36238444 24 HOM del
- 16642 chr6 36238445 24 HOM del
- 16642 chr6 36238446 24 HET del
- 16642 chr6 36238447 24 HET del
- 16642 chrX 2336324 24 HET del
- 16642 chrX 2336325 24 HOM del
- 16642 chrX 2336326 24 HOM del
- 16642 chrX 2336327 24 HOM del
- [root@localhost ~]# ./p | column -t
- 16642 chr6 36238437-36238445 HOM del NM_001145717 PNPLA1 CDS1:intron1
- 16642 chrX 2336325-2336326 HOM del NM_002364 MAGEB2 5-UTR
- [root@localhost ~]#
复制代码 |
|