1. Convert your SAM file to BAM file.
2. Sort the BAM file using chromosomal coordinates.
3. Then try: samtools view -F 4 Input.bam | cut -f3 | uniq -c

You can use samtools (http://samtools.sourceforge.net/samtools.shtml)

Or you can try awk command on your SAM file:

awk '{if ($5 > 1)cnt[$3]++}END{for (x in cnt){print x,cnt[x]}}'  Input.sam

to count unique read names per contig, you can first extract unique pairs <read> <contig> and then count only contig names:

samtools view file.bam | cut -f1,3 | sort | uniq | cut -f2 | sort | uniq -c

This will yield a list of pairs <number of reads> <contig>

#!usr/bin/perl-w
use strict;
use warnings;

my %hash = ();
while(<>){
        chomp;
    next if($_ =~ /^@/); ## remove the headers in sam file
        my @s =split;
        my $contig = $s[2];
        $hash{$contig}++;
}

foreach (sort keys %hash){
        print join ("\t", $_,$hash{$_}),"\n";
}

Run the above code on the sam file as

perl read_count.pl test.sam > counts.txt

Hope this helps

Varun

I might understand this wrongly but if you have a BAM file (which was indexed with samtools index) you can just run

samtools idxstats your.bam

which will list the "name, sequence length, # mapped reads and # unmapped reads" for each sequence/chromosome (see the idxstats section at http://samtools.sourceforge.net/samtools.shtml)

Andreas

Edit: Missed the fact that you were looking for unique counts

Hi Varun,

Thanks a lot.

Yes It is working,but somehow it is counting 1 time less than the actual number of times it is mapped.Just to be clear if a contig is mapped to "29" ID's it counts only "28" and is it possible to print the "ID's" that were mapped along as the 3rd column with the total number like now.

Hi Pavel Senin ,

Thank you for the reply.

Actually what you meant in the first post was right. Regardless of how many times a read is mapped through out the file, at the end it should print "contig name" and total number of reads mapped to that particular contig along with the list of read ID's.

Hi ashutoshmits,

Thank you for the answer.

I didn't try with the first option yet,but with "awk" it doesn't gave the total count correctly.