intersect of two sam files

[go2seq]

I used two aligners for paired end reads, and got two sam files. Are there any ware or scripts to implement these two operations?

- How to find the intersect of two sam files by read ids? To use bedtools, one needs to convert sam to bam to bed and work on 2 bed files.

- I'd like to merge the 2 sam files and, for duplicated reads being aligned in both, select a better alignment by the number of mismatch (NM).

Please advise.

[ShaunMahony]

samtools has a merge function, but I don't know how well it performs, and the files have to be sorted before merging. You could always use samtools to convert to SAM (if you have BAM), concatenate the files from the command line (making sure that the header only goes once to the top of the file), and using samtools again to sort.

[swbarnes2]

Actually, bedtools will take bam files for some applications.

[malcook]

Here is a strategy, assuming that

each sam file has at most one alignment per read
the sam files are not (yet) sorted

Use `samtools view` to extract just the alignments from each file

Use perl (or something) the extract the value of NM and pre-pend it to each line as a new first column.

`sort` the two alignments together first by readname (-k2,2) and secondarily descending by score (-k1,1nr)

`cut` the NM (1st column) score out of your sorted/combined data

use `sort` again to pick out the first (best scoring) alignment in each group (by using the --merge --unique -k1,1)

somehow contrive to reheader this combined set of alignments