ha thanks Brent, I didn't know this syntax for awk.
• 0 views
•
link
written
by
Pierre Lindenbaum
What Is The Most Direct Way To Extract The Splice Junctions From A Sam File?
Hi!
I need to get the introns coordinates (chr:strat-end and the strand) of the spliced reads wich are in SAM file.
I have no experience with this kind of format, so I plan to transform the SAMs files to BED12 files with this pipe line:
samtools view -bS -o file.bam file.sam && bamToBed -bed12 -i file.bam > file.bed12
And then, use galaxy or own python script (which I haven't wrote yet) to extract the introns from the spliced reads.
The problem with this method is weight of the files... it'll take some time to convert to BED12 and then extract the introns...
Do you know a more direct way to extract the introns coordinates from the SAMs files?
2 answers
It's not clear exactly what you want to do, but you can extract only the spliced alignments with:
samtools view -S file.sam | awk '($6 ~ /N/)'
and you may be able to pipe that command directly to bamToBed as:
samtools view -hS file.sam | awk '($6 ~ /N/)' \
| samtools view -bS - \
| bamToBed -bed12 > file.bed12
From there, you can get introns by subtracting exons:
so:
bed12ToBed6 -i file.bed12 > file.bed6
subtractBed -a file.bed12 -b file.bed6 -s | cut -f 1-6 > file.introns.bed
Hi brentp!... I not only want to get the spliced reads... I want to get the genomics coordinates of the introns in those reads, I mean where they are in the genome (chromosome, start, end and strand). From a bed12 you can easily extract this information, but I wonder if are a script to extract direct from the SAM file.
if you swing with R, a slight change to my post to account for strandedness should get you there...
• 0 views
•
link
written
by
Malcolm.Cook
Log in to answer this question.
More posts like this
-
Bamtools : Extract Sense And Anti-Sense Reads
written by Nicolas Rosewick
-
No Coverage Information In *.Fpkm_Tracking Files Obtained From Cuffdiff Output
written by komal.rathi
-
Calculating Coordinates Of Split Rna-Seq Reads Converted From Bam To Bed
written by Alex Reynolds
-
Samtools Pipelines
written by Geparada
-
Input Error Using Varscan2.3
written by Jordan
-
Emule "Samtools View" Using Pysam
written by Geparada
-
Fastq With Format Errors
written by Geparada
-
Liftover Of Na12878 Coordinates Into Hg19 Coordinates
written by Geparada
-
Gencode Gtf To Bed12 Or Psl
written by Geparada
-
Aligning Colorspace Reads Using Bwa
written by Farhat
The intron coordinates aren't really there, if anything, they are implicit in (many of) the alignments. Are you sure you want to roll your own instead of running tophat or mgene or similar?
Yes, I know there are implicit in bed12... But I'm not worry about that, because I wrote a scrip to extract introns coordinates from psl alignment format, so I'll try to adapt to bed12 inputs, or do as brentp proposes...