4
$\begingroup$

I am interested in identifying indels in whole genome bisulfite sequencing data (76bp paired end). Currently, I do this by setting the -rfg and -rdg affine gap penalty scores for bowtie2 to more permissive values than the default 5+3N and mapping using the bisulfite sequencing alignment wrapper, bismark.

My question is, what values of -rfg and -rdg will allow me to identify longest possible indels without sacrificing alignment quality? Is it better to set the affine penalty to zero with a high penalty for initially opening an indel (ex. 8+0N)? Or is it better to keep the initial penalty low and having a nonzero penalty for extension (ex. 2+1N)?

Improve this question
$\endgroup$

1 Answer 1

Reset to default
3
$\begingroup$

You'd be best off by starting with -rfg and -rdg as is and reverting bismark's change of --score-min back to the default for bowtie2. That alone will allow for much longer indels. If that still doesn't suffice, then I'd play around more with --score-min before messing with the gap open/extend penalties. If you do need to play with those, then increase the gap open penalty and decrease the gap extension to 1. Do check that the resulting alignments aren't nonsense though!

Improve this answer
$\endgroup$
1
  • $\begingroup$ With the bismark default --score-min, -rdg, and -rfg settings, it appears that the maximum sized indel identifiable with my read length is 3bp. I changed both -rdg and -rfg to a gap open/extension of (1,1) and the maximum length I can detect is 14bp. I will try modifying --score-min only. It's a fine line to walk, detecting the longest possible insertions while keeping the alignment clean! $\endgroup$
    – Ben D.
    Commented Aug 4, 2017 at 3:38

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.