This is a test version of Biostars. For the public version, visit https://www.biostars.org.
Removing chrEBV reads from bam file

Hello!

I am trying to visualize on the UCSC Genome Browser some .wig files, but I keep getting the following error message during the uploading of my custom track:

Error line 39996036 of somefile.gz: 'chrEBV' is not a valid sequence name in hg38

I have already trying to filter out from the starting .bam files all the reads corresponding to chrEBV, using the procedure suggested in Remove mitochondrial reads from BAM files.

Unfortunately, I still get the same error when trying to upload the new .wig file.

How can I make sure to get rid of all the chrEBV reads and finally manage to visualize my data on the Genome Browser?

Thanks!

bam genome
• 4,227 views • link
updated by Ram Ram • written by Marco Pannone • Marco Pannone

I have already trying to filter out from the starting .bam files all the read

how did you do that ? what was the cmd -line ? what is the output of `samtools idxstat' ?

• 0 views • link
written by Pierre Lindenbaum Pierre Lindenbaum

My command line was: samtools idxstats input.bam | cut -f 1 | grep -v chrEBV | xargs samtools view -b input.bam > output_filtered.bam

• 0 views • link
written by Marco Pannone • Marco Pannone

I first created .bigWig files with the following command line:

bamCoverage -b input.bam -o output.bw -of bigwig -bs 20 -p 6 --effectiveGenomeSize 2747877777 --normalizeUsing RPKM -e 76 --centerReads

and then converted it from .bigWig to .wig, using this command line:

bigWigToWig input.bw output.wig

bamCoverage does not allow me to create a .wig file directly, that's why I follow this two-steps procedure. Also, when uploading the tracks on the UCSC Genome Browser, I prefer to go with .wig.gz files rather than .bigWig, since for .bigWig I need to upload them first on a web-server and then provide an URL (as long as I understood).

• 0 views • link
written by Marco Pannone • Marco Pannone

yes, I think ATPoint is right. If bamCoverage uses the SAM header dictionary to create the wig, i will create some records containing the chrEBV chromosome.

• 0 views • link
written by Pierre Lindenbaum Pierre Lindenbaum

Yes, indeed his approach worked fine. Thanks anyway for the reply and the time dedicated to my issue :)

• 0 views • link
written by Marco Pannone • Marco Pannone

1 answer

That will only remove the sequences aligned to chrEBV but not remove it from the header. Depending which tool you use there will be chrEBV in the wig with coverage of zero.

Try (untested):

samtools idxstats input.bam \
| cut -f 1 \
| grep -v 'chrEBV' \
| xargs samtools view -h input.bam \
| grep -v 'chrEBV' \
| samtools view -o output_filtered.bam
• 0 views • link
written by ATpoint ATpoint

This worked for me! Now the .wig files upload on the UCSC Genome Browser without any error message. Thanks a lot for that!

Just out of curiosity, how is it possible that I ended up having chrEBV reads in my dataset? They belong to the Epstein-Barr virus, right? The biological starting material of the experiment is human cell culture, so I wonder if this is a sign of cell culture contamination. Still, I did not expect that GRChM38 genome assembly would include also genomic regions belonging to EBV.

• 0 views • link
written by Marco Pannone • Marco Pannone

Adding on this, here is a nice explanation as well: C: What is the story behind chrEBV?

It basically captures EBV contaminations avoiding that these sequence false-positively get aligned to the genome by providing the true-positive sequence (=the chrEBV). This is a so-called "decoy" sequence.

• 0 views • link
written by ATpoint ATpoint

Log in to answer this question.

More posts like this