Hello,
I just aligned my first piece of RNA-Seq data and have a .bam file. I'd like to make a file of counts corresponding to genes that I can use in DESeq2. I downloaded a list of all genes from the UCSC gene browser and saved it as a .bam file. I'm trying to do the following:
bedtools multicov -bams UCSC_knownGene.bam -bed myData_out/accepted_hits.bam
But then I get this error:
Could not open input BAM files.
I've only read about this error on forums where people are using multiple bam files, so I'm not sure why I'm getting it if I'm only submitting one file that was directly downloaded from UCSC Genome Browser. Has anyone had a similar issue? Is there a better way to get .bam files to use for this?
Thank you!
•
written
by
Kristin Muench
1 answer
You immediate problem is likely solved by:
bedtools multicov -bams myData_out/accepted_hits.bam -bed UCSC_knownGene.bam
I'm assuming you just misnamed the knownGene file and that it's really a bed12 file.
Your mediate problem is that you shouldn't do any of this to begin with. Bedtools multicov should not be used to get counts for RNAseq. The resulting counts violate the assumptions made by tools like DESeq2 (and edgeR and limma and many more). Use featureCounts or htseq-count and download a GTF file from UCSC.
Log in to answer this question.
-
Visualizing alignment of kallisto .bam files to hg19 (and custom plasmids)
written by Kristin Muench
-
bedtools multicov displays only the last BED feature
written by ThePresident
-
how to calculate coverage in bam file over a specific bed file
written by Lila M
-
List of Ensembl Transcript IDs corresponding to protein-coding genes
written by Kristin Muench
-
Using Tophat2 to align illumina files to ucsc hg19
written by fk566938 •
-
[Solved] Problems merging TopHat mapped and unmapped BAMs for submission to ENA
written by Chris Cole
-
bedtools multicov and fraction overlap flag
written by fitziano •
-
Extracting a trimmed output from a bam file
written by Will
-
Using Gnu Parallel For Bedtools
written by GouthamAtla
-
Bedtools Multicov Need A Bam Index File Specification Option
written by ewre
Your list of genes should be in bed format.
https://genome.ucsc.edu/FAQ/FAQformat.html#format1
Also is your bam file indexed? Use
samtools indexto index your bam fileNot sure it was indexed - thank you!
Try:
it will probably say "ASCII text". If it a text file, then do:
If unsure, post here the result of head. Bam files should be gzipped. I suspect the file you downloaded gene list is either bed (as Varun Gupta pointed out) or gtf format.