This is a test version of Biostars. For the public version, visit https://www.biostars.org.
How to obtain insert sequence from paired-end mapping result (SAM/BAM)

I have paired-end reads R1 and R2.

|-----R1-----|----------------------------|-----R2-----|
       A                     B                  C

Total insert here is A + B + C. How can I get this sequence from the BAM/SAM output?

The reads are mapped using classical Bowtie2 parameter.

    bowtie2 -x myindex \
   -1 R1.fastq \
   -2 R2.fastq \
   --no-discordant \
   --no-mixed -p 4 \
    --very-sensitive  \
   --no-mixed \
    | samtools view -bS - > output bam
rna-seq bowtie2 alignment samtools
• 3,530 views • link
updated by Fabio Marroni Fabio Marroni • written by gundalav gundalav

You can't really get that solely from a SAM/BAM file, since the sequence appropriate for B for a particular fragment is unknown. What is your actual biological goal?

• 1 views • link
written by Devon Ryan Devon Ryan

If we know the mapping location of A and C, we can get A+B+C there right? With Picard we can get the A+B+C length using pathtopicard\CollectInsertSizeMetrics.jar I wonder if there is any existing tool to actually extract the fragment from the genome. At the end of the day want to translate A+B+C into protein and derive MHClass-I core from there.

• 1 views • link
written by gundalav gundalav

It's unclear if you want the sequence of the fragment or the sequence of the reference genome covered by the fragment. For the latter, sure, that's easy enough with a small bit of scripting (I'd use pysam in python, but to each his/her own). For the former, see my earlier comment.

• 1 views • link
written by Devon Ryan Devon Ryan

1 answer

Ok, if your insert strt at ChrXXX position Astart and end on the same ChrXXX at position Cend.

Then you create a bedfile and you call it yourbed.bed

ChrXXX Astart Cend

and then use

bedtools getfasta -fi yourfasta.fa -bed yourbed.bed

Done!

However, I re-emphasize what Devon just said: you do not know what the real sequence of B is. You are assuming that your sample has identical sequence to the reference.

• 1 views • link
written by Fabio Marroni Fabio Marroni

Log in to answer this question.

More posts like this