Filter bam file based on cigar string?

Question

[bioinf_sci]

How can I filter a bam file using samtools view that retains reads where the cigar string contains M<=35?

Answer 1

[jkbonfield]

A bit messy, but it can be done with regexp matching:

samtools view -e 'cigar =~ "(^|[^0-9])([0-9]|[12][0-9]|3[0-5])M"' -o out.bam in.bam

The regexp is start or non-digit, followed by ?M [12]?M or 3[012345]M.

[Ruoshui]

Hi, if I want to discard reads where the cigar string contains M < 100, how should I write this instead? Thank you~

Answer 2

[Pierre Lindenbaum]

using samjdk: https://jvarkit.readthedocs.io/en/latest/SamJdk/

I cannot post code here ! there is a bug ni biostars...

so here is a gist:

java -jar src/jvarkit-git/dist/jvarkit.jar samjdk -e 'return !record.getReadUnmappedFlag() && record.getCigar().getCigarElements().stream().filter(CE->CE.getOperator().equals(CigarOperator.M)).anyMatch(CE->CE.getLength()<=45);' in.bam

view raw README.md hosted with ❤ by GitHub

[GenoMax]

java -jar src/jvarkit-git/dist/jvarkit.jar samjdk -e 'return !record.getReadUnmappedFlag() && record.getCigar().getCigarElements().stream().filter(CE->CE.getOperator().equals(CigarOperator.M)).anyMatch(CE->CE.getLength()<=45);' in.bam

Posting code here in case gist goes away.