Removing contaminating reads

[travelk]

Hey all,

I have a very basic question that I just cannot seem to find a straight answer to despite scouring Google and SEQanswers.

When I ran my dataset through Fastq_screen I discovered I have E. coli contamination. I ran my files in Tophat2 against the E.coli genome and determined that about 5% of my reads positively map to E. coli.

I would like to remove these reads but I cannot figure out how.

Most forum posts suggest to remove them but do not specify the exact method to do that. My initial thought was to take my unmapped.bam file from my E.coli alignment and use that as my "clean" data but it means converting the data back to fastq in paired reads format and adding various processing steps so I'm not really sure what is left in that file. Naturally, I'd like to lose as little data as possible.

Is there a simple way to just take my original file fastq file, extract the E.coli+ reads and have a clean data set?

Or alternatively, does having a contamination like that not really matter because they won't map to the mouse genome anyway?

Obviously I'm pretty new to this and teaching myself so any advice, however simple it may seem, is much appreciated.

Thanks for your help!

[dpryan]

Just align to the E coli genome with bowtie2, which can be told to write unmapped reads in fastq format. You can then directly align the resulting fastq files.

[a.kmg]

You can run Bowtie2 with your fastq file on E.coli to suppress reads corresponding to E.coli :

bowtie2 -U rawFastq -p nbproc --un fasqtFileWithoutEcoli indexEcoli -S ecoli_reads.sam

This command treats the raw fastq file, creating a new fastq containing not E.coli reads (--un = unmapped reads) and recovering reads aligned on E. coli in a sam file (-S option). Delete -S option if you do not want to get the E.coli reads.

[GenoMax]

BBMap example (remain in the same directory or adjust paths accordingly). Will work on PC/Mac/*nix.

1. Get Ecoli genome fasta file.

2. Build an index for Ecoli genome using BBMap.

$ bbmap.sh ref=./E_coli_genome.fa

3. Align against the Ecoli genome index saving reads that don't align to a new file (these are the reads you want, the outu file below)

$ bbmap.sh in=your_fastq_file path=./ outm=E_coli_reads.fastq outu=reads_you_want.fastq qin=33

[travelk]

Thank you! This totally did the trick. I processed the files and then reran them through fastq_screen and all the E.coli reads were gone.

And thank you a.kmg and GenoMax for your clear step by step instructions.

I ended up using bowtie2 because that it what I've been working with so far. I initially got an error message but realized I was missing the -x in the command line. So, in the end I used:

bowtie2 -U raw.fastq -p nbproc --un FileWithoutEcoli.fastq -x indexEcoli -S ecoli_reads.sam

[nareshvasani]

Hi fellow,

I am trying to find rRNA gene from my input file:

Step1: Create Index

#bowtie2-build rRNA.fasta rRNA.index

Step 2: Align to rRNA index inorder to get rRNA free fastq file.

#bowtie2-align -p 2 -k 1 -q -U /filter_clean.fastq --un fasqFileWithoutrRNA -x rRNA.index

When I run second step, it comes with error saying:

" bowtie2-align: option '--un' is ambiguous; possibilities: '--ungapped' '--unpaired' "

So I replace --un with --unpaired, But it is not working.

Can anyone please shed some light on this.

I would really appreciate your help.

Thanks,
Naresh

[jpnm]

As far as I know "--un" is not an option of bowtie2-align...but only of bowtie2!

So if you run, it should work fine!

bowtie2 -p 2 -k 1 -q -U /filter_clean.fastq --un fasqFileWithoutrRNA -x rRNA.index

I hope it helps!

[nareshvasani]

As far as I know "--un" is not an option of bowtie2-align...but only of bowtie2!

So if you run, it should work fine!

bowtie2 -p 2 -k 1 -q -U /filter_clean.fastq --un fasqFileWithoutrRNA -x rRNA.index

I hope it helps!

Hi jpnm,

Thanks for your input, I am trying you suggestion. Will keep you posted.


Naresh