Hi guys,

I need to identify insertions of specific sequences in my genomes. I have WGS data with read length of 250 nucl. Previously I was working on identification of large sequences insertions and simply looked for discordantly mapped reads for that using bowtie2 and bbmap. 
Now I am looking for sequences that are shorter than read length. As I want these data to be comparable to my previous computations I want to use bowtie2 for alignment again. As far as I understand, I should use local mode. I've tried this, but got alignment of very short parts of my reads. Is there a possibility to require the alignment length to be equal to reference sequence length? When I identify the reads that are mapped to the reference sequences, I will find where are their pairs mapped to the reference genome to find the insertion site, is it a correct strategy?

Thanks id advance!
You can call insertions near read length with BBMap like this (assuming paired interleaved reads that mostly overlap):

    bbmerge.sh in=reads.fq out=merged.fq
    reformat.sh in=merged.fq out=filtered.fq minlen=350
    bbmap.sh in=filtered.fq out=mapped.sam ref=ref.fa slow
    callvariants.sh in=mapped.sam out=vars.vcf ref=ref.fa ploidy=2

That's a variant of a longer and slower pipeline that I used to successfully call insertions from 1bp to 219bp using 2x100bp reads, on E.coli.  If you have sufficient time, coverage, and memory, you can use this method (adjusted for 2x250bp reads) instead for even longer insertions (it's great on small genomes but would take a lot of memory for human-sized ones):

    filterbytile.sh in=reads.fq.gz out=filtered_by_tile.fq.gz
    bbduk.sh in=filtered_by_tile.fq.gz out=trimmed.fq.gz maxns=0 ktrim=r k=23 mink=11 hdist=1 tbo tpe minlen=220 ref=adapters.fa
    clumpify.sh in=trimmed.fq.gz out=eccc.fq.gz passes=4 ecc unpair repair
    tadpole.sh in=eccc.fq.gz out=ecct.fq.gz ecc prefilter=2
    tadpole.sh in=ecct.fq.gz out=extended.fq.gz ordered mode=extend el=20 er=20 k=62 prefilter=2
    bbmerge-auto.sh in=extended.fq.gz out=merged.fq.gz outu=unmerged.fq.gz rem k=93 extend2=120 prefilter=2
    bbmap.sh in=merged.fq.gz out=merged.sam.gz slow ref=ref.fa
    callvariants.sh in=mapped.sam out=vars.vcf ref=ref.fa ploidy=2

