SEQanswers

Go Back   SEQanswers > Bioinformatics > Bioinformatics



Similar Threads
Thread Thread Starter Forum Replies Last Post
intersect (actually: filter) a gtf file with coordinates from a bed-file dietmar13 Bioinformatics 6 07-13-2017 06:20 AM
Sorted BAM read count >2x total FASTQ count pkMyt1 Bioinformatics 8 05-13-2015 10:43 AM
How to extract read count from bam file lylyf1987 RNA Sequencing 0 04-07-2015 12:19 PM
Getting mapped read count from a BED file - bedtools coverage? samhokin RNA Sequencing 0 12-20-2013 02:22 PM
BEDtools intersect output is BED instead of BAM syfo Bioinformatics 1 12-18-2012 05:26 AM

Reply
 
Thread Tools
Old 06-13-2016, 02:42 AM   #1
edeka
Member
 
Location: us

Join Date: Apr 2016
Posts: 14
Default intersect file bam-bed and count read

Hi all!
I have RNA-seq data and I used STAR in order to align the reads.
I would like to use my bed file in order to have the only read that are mapped in the position reported in bed file and after I would like to know how many read are mapped in each interval.
For the first part I tried to use "bedtools intersect":
bedtools intersect -abam -a file.bam -b file.bed
I saw that the number of the read reported in the output file were less than the input file but It doesn't seem right. In fact, for example, in my bed file for chromosome the first interval begins with the position 335 and in my new output file there are reads mapped in position 80.

Thank you in advance and I am sorry for my probably bad English!
Best
edeka is offline   Reply With Quote
Old 06-13-2016, 06:08 AM   #2
GenoMax
Senior Member
 
Location: East Coast USA

Join Date: Feb 2008
Posts: 7,087
Default

The option you should look at is coverageBed.

As for the other observation if you are providing a set of BED intervals then only reads that fall in those intervals will be reported (so that number can be expected to be less than the input). It is possible that the read in question is mapped beginning at position 80 but depending on how long it is it must be extending into the first interval.
GenoMax is offline   Reply With Quote
Old 06-14-2016, 12:43 AM   #3
SylvainL
Senior Member
 
Location: Geneva

Join Date: Feb 2012
Posts: 179
Default

Since you did the mapping with STAR, it is possible a read will be really long once mapped to the genome if it is at an exon-exon junction... You can use the parameter -split in your bedtools command to report only the region your read is really mapping...
SylvainL is offline   Reply With Quote
Old 06-20-2016, 03:34 AM   #4
edeka
Member
 
Location: us

Join Date: Apr 2016
Posts: 14
Default

Quote:
Originally Posted by GenoMax View Post
The option you should look at is coverageBed.

As for the other observation if you are providing a set of BED intervals then only reads that fall in those intervals will be reported (so that number can be expected to be less than the input). It is possible that the read in question is mapped beginning at position 80 but depending on how long it is it must be extending into the first interval.
Thank you for your reply. The read is about 20 nucleotides, so beedtools should discard this read...
edeka is offline   Reply With Quote
Reply

Tags
bedtools intersect, count, star aligner

Thread Tools

Posting Rules
You may not post new threads
You may not post replies
You may not post attachments
You may not edit your posts

BB code is On
Smilies are On
[IMG] code is On
HTML code is Off




All times are GMT -8. The time now is 08:46 AM.


Powered by vBulletin® Version 3.8.9
Copyright ©2000 - 2021, vBulletin Solutions, Inc.
Single Sign On provided by vBSSO