Hi biostars,
By tophat2 (tophat2 -p 8 -G genes.gtf genome file.fastq) command (gtf I think is the annotation file and genome is the whole genome fasta), I produced a file named accepted_hits.bam which using
bam2bed < accepted_hits.bam | bedmap --echo --count genes.bed - > answer4.bed
command, by bedops tool or bedtools multicov -bams yourBAMfile.bam -bed yourGenes.bed command
by bedtools, I have a bed file now, like below
I 334 646 "YAL069W . + protein_coding CDS 0 exon_number "1"; gene_id "YAL069W"; gene_name "YAL069W"; p_id "P3633"; protein_id "YAL069W"; transcript_id "YAL069W"; transcript_name "YAL069W"; tss_id "TSS1128"; 2**
For example 2 reads have been mapped on gene YAL069W in 334-646 distance
but how I can know where each of these reads individually have been mapped? I mean how I can have a file in which I could see where each read individually has been mapped. For example these two reads in protein_coding (CDS) part, where have been mapped individually not just knowing that in 334 - 646, two reads have been mapped on gene YAL069W (not a distance in where a number of reads have been mapped).
Please if you have any idea, let me know
Thank you so much
1 answer
You can view your bam files in IGV to see to positions of where each read have mapped.
It only requires indexing of your bam file, e.g. using samtools:
samtools index accepted_hits.bam
•
written
by
Annika Forsingdal
Log in to answer this question.
-
Putting bedtools in the path
written by zizigolu
-
Distribution of mapped reads
written by ilovesuperheroes1993 •
-
How to get regions with mapped reads
written by nanoide
-
Extracting mapped paired end FastQ reads from BAM file using Samtools - trying to understand -f 0x0…
written by acxcrew007 •
-
Problem converting bam to sam (samtools)
written by hinkel2 •
-
Tophat unable to produce accepted_hits.bam file for single end read data
written by nikhilvgbt •
-
the times each genes have been repeated in my file is equal to the number of reads have been mapped…
written by zizigolu
-
the number of reads that have mapped to the each gene
written by zizigolu
-
why the bed produced by bowtie2 and tophat2 are different?
written by zizigolu
-
How To Handle 'Mate Mapped To Different Chr' In Rna-Seq Data Analysis ?
written by GouthamAtla
You never encountered SAM format in your life? Samtools does not ring the bell? RTFM
Hey don't be angry please. If you read my previous posts you could observe that I tried bedops, bedtools, samtools. I have a file with number of reads but this number is the sum of mapped reads in a given distance but I need to know where the each read is mapped individualy. I used samtools
For example
samtools view yourBAMfile.bam chr1:567876-568100 | wc -lbut there is no output just running something in cmd. After running this is the example of output of this command in cmdbut again this is just the total number of the mapped reads in given distance and not clear the exact position of where each read has been mapped
If you expect that piping
samtools viewoutput intowcwill give you where the individual reads are mapping then you need to go back and learn what these commands do. Including being able to check that you can execute samtools, that the BAM file is there, it is not corrupted/overwritten by some badly executed command etc. Nobody can help you with this except by dictating line by line and resolving problems people have in their first few hrs of working with a command line. Read e.g.:http://openwetware.org/wiki/Wikiomics:WinterSchool_day1#Introduction_to_Linux_and_the_command_line
And check that you have samtools on the PATH & chr1 in your BAM file header SQ as a sanity check.
Good luck
Since I joined biostars never I have faced such an angry biostar! I'm working with linux since about a month ago. Anyway thanks
Hello Fereshteh!
We believe that this post does not fit the main topic of this site.
Not a slightest attempt at solving the problem before posting
For this reason we have closed your question. This allows us to keep the site focused on the topics that the community can help with.
If you disagree please tell us why in a reply below, we'll be happy to talk about it.
Cheers!