Compare commits

..

No commits in common. "4443b3c2a939b25d0eebbe8b2d1003a47501bc3d" and "55658e204f8b01710d96b670f2b3eeb3274193d4" have entirely different histories.

View File

@ -2,7 +2,6 @@
use strict; use strict;
use Parallel::ForkManager; use Parallel::ForkManager;
# this script creats one file with UMI unique reads and one with UMI duplicated reads # this script creats one file with UMI unique reads and one with UMI duplicated reads
# as input you need aligned sorted by coordinate bam file
my @torun = (); my @torun = ();
foreach my $file ( <*Aligned.sortedByCoord.out.bam> ) { foreach my $file ( <*Aligned.sortedByCoord.out.bam> ) {
push @torun, $file; push @torun, $file;
@ -34,12 +33,12 @@ foreach my $file ( @torun ) {
my $uniq = join( ':', $chr, $pos, $flag, $tlen, $bc ); my $uniq = join( ':', $chr, $pos, $flag, $tlen, $bc );
my $pos_ = $pos-1; my $pos_ = $pos-1;
while ( $cigar =~ m/(\d+)([SHMDIN=])/g ) { while ( $cigar =~ m/(\d+)([SHMDIN=])/g ) {
$pos_+=$1 if $2 eq 'M' or $2 eq '=' or $2 eq 'D' or $2 eq 'N'; # find position for minus strand $pos_+=$1 if $2 eq 'M' or $2 eq '=' or $2 eq 'D' or $2 eq 'N';
} }
my $uniq2 = join( ':', $chr, $pos_, $flag, $tlen, $bc ); my $uniq2 = join( ':', $chr, $pos_, $flag, $tlen, $bc );
if ( exists($duplicates{$id}) or # already marked as duplicate if ( exists($duplicates{$id}) or # already marked as duplicate
( not($flag&16) and ++$seen{ $uniq } > 1 ) or # plus strand ( not($flag&16) and ++$seen{ $uniq } > 1 ) or # plus strand
( $flag&16 and ++$seen{ $uniq2 } > 1 ) #minus strand ( $flag&16 and ++$seen{ $uniq2 } > 1 )
) { ) {
print F2; print F2;
$dups++; $dups++;