Answers.Body: <p>I was asked to write some code sometime near the end of last year to do this, so that viral sequences could be split based on variants at a particular location. The code I wrote doesn't actually do the splitting - it creates SAM read groups based on those variants - but that splitting can be done after processing with this script using <code>samtools split</code>.</p>
  <p>As with most of my scripts, it hasn't had much testing / use:</p>
  <p><a href="https://gitlab.com/gringer/bioinfscripts/-/blob/master/samVarSplitter.pl" rel="nofollow noreferrer">https://gitlab.com/gringer/bioinfscripts/-/blob/master/samVarSplitter.pl</a></p>
  <p>usage:</p>
  <pre><code>samtools view -h mapped_reads.bam | ./samVarSplitter.pl [-ref &lt;ref&gt;] [-pos &lt;int&gt;] [options]
  </code></pre>
  <p>example:</p>
  <pre><code>samtools view -h mapped_reads.bam | ./samVarSplitter.pl -ref Wuhan-Hu-1 -pos 11083
  </code></pre>
  <p>If the input file contains a lot of reads that don't include the variant location, I'd recommend pre-filtering reads to only include those that cover that location, because expanding the CIGAR string and storing a record of the positions can take some time.</p>
  <p>Note: currently the position argument only supports a single location. I could extend this to multiple positions (i.e. a position range) if that were a desired feature.</p>
  
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2020-08-31T18:54:29.273Z
Answers.Text: Thanks so much, this is pretty much what I was looking for! If I understand it correctly, the script only works with a single position though correct? The deletion I'm looking at is 9bp in size, but I suppose I can just use one of the positions as a stand-in for the entire deletion.
Answers.UserDisplayName: jazzbo
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2020-08-29T20:58:37.440Z
Answers.LastActivityDate: 2020-09-02T04:05:31.463Z
Answers.LastEditDate: 2020-09-02T04:05:31.463Z
Answers.UserDisplayName: gringer
Answers.Body: <p>The CIGAR string should tell you if a read has a deletion in it.  So filter the bam files based on that.</p>
  
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2020-08-27T08:26:47.940Z
Answers.Text: The OP is looking for a known deletion, so if you look to identify reads in the region with 'D' in the cigar, you might include other deletions
Answers.UserDisplayName: StupidWolf
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2020-08-27T00:43:10.707Z
Answers.LastActivityDate: 2020-08-27T00:43:10.707Z
Answers.UserDisplayName: swbarnes2
Body: <p>I have an indexed BAM file containing long-read sequencing data and I'd like to split the reads contained within into those with a known deletion and those without the deletion (I have the deletion coordinates available to me) when mapped against the hg38 genome. Ideally I'd like to end up with two BAM files (one with reads containing the deletion and one without) for downstream analysis. I've searched a bit online and haven't really found any tools that do exactly that, but this seems like a common enough task that there has to be something.</p>
  <p>Thanks in advance for any help/guidance!</p>
  <p><strong>UPDATE</strong></p>
  <p>In addition to @gringer's great answer, I decided to write a small python script myself that solves my issue using pysam, the fruits of which can be found <a href="https://github.com/brookshenry3/bam_splitter" rel="nofollow noreferrer">here</a>. If anyone has any feedback or suggestions on it I'll happily take them, as I'm still very new to this. It's undergone very little testing, so user be warned.</p>
  
Tags: <sam><long-reads><indel>
Title: Sort reads in BAM file based on presence of specific deletion?
ThreadText: # Sort reads in BAM file based on presence of specific deletion?
  
  (Asked by: username_0 on 2020-08-26)
  
  I have an indexed BAM file containing long-read sequencing data and I'd like to split the reads contained within into those with a known deletion and those without the deletion (I have the deletion coordinates available to me) when mapped against the hg38 genome. Ideally I'd like to end up with two BAM files (one with reads containing the deletion and one without) for downstream analysis. I've searched a bit online and haven't really found any tools that do exactly that, but this seems like a common enough task that there has to be something.
  
  Thanks in advance for any help/guidance!
  
  **UPDATE**
  
  In addition to @username_1's great answer, I decided to write a small python script myself that solves my issue using pysam, the fruits of which can be found [here](https://github.com/brookshenry3/bam_splitter). If anyone has any feedback or suggestions on it I'll happily take them, as I'm still very new to this. It's undergone very little testing, so user be warned.
  
  **Tags:** sam, long-reads, indel
  
  
  ---
  
  ## 2 answers
  
  ### Answer by username_1 on 2020-08-29. Score: 1 | Accepted answer
  
  I was asked to write some code sometime near the end of last year to do this, so that viral sequences could be split based on variants at a particular location. The code I wrote doesn't actually do the splitting - it creates SAM read groups based on those variants - but that splitting can be done after processing with this script using `samtools split`.
  
  As with most of my scripts, it hasn't had much testing / use:
  
  <https://gitlab.com/username_1/bioinfscripts/-/blob/master/samVarSplitter.pl>
  
  usage:
  
  ```
  samtools view -h mapped_reads.bam | ./samVarSplitter.pl [-ref <ref>] [-pos <int>] [options]
  
  ```
  
  example:
  
  ```
  samtools view -h mapped_reads.bam | ./samVarSplitter.pl -ref Wuhan-Hu-1 -pos 11083
  
  ```
  
  If the input file contains a lot of reads that don't include the variant location, I'd recommend pre-filtering reads to only include those that cover that location, because expanding the CIGAR string and storing a record of the positions can take some time.
  
  Note: currently the position argument only supports a single location. I could extend this to multiple positions (i.e. a position range) if that were a desired feature.
  
  > **Comments:**
  > **username_0**: Thanks so much, this is pretty much what I was looking for! If I understand it correctly, the script only works with a single position though correct? The deletion I'm looking at is 9bp in size, but I suppose I can just use one of the positions as a stand-in for the entire deletion.  
  
  ---
  
  ### Answer by username_2 on 2020-08-27. Score: 2
  
  The CIGAR string should tell you if a read has a deletion in it. So filter the bam files based on that.
  
  > **Comments:**
  > **username_3**: The OP is looking for a known deletion, so if you look to identify reads in the region with 'D' in the cigar, you might include other deletions  
  
  ---
  
  