Hi,
As part of my chip seq analysis, I tried to run a script to convert fastq file into .bam file without the creation of a .sam file (using piping). This is the script:
${bowtie2_source} -x ${ref_genome} -U ${fastq_file} -S | ${samtools} view -bS - ${target_dir}/${sample_name}.bam
unfortunately, I recieved the following error: [main_samview] random alignment retrieval only works for indexed BAM or CRAM files. How can I fix this?
•
written
by
Bella_p
you wrote:
${samtools} view -bS - ${target_dir}/${sample_name}.bam
you want
${samtools} view -bS -o ${target_dir}/${sample_name}.bam -
Great, that works perfectly. Excuse me but I'm a begginer and trying to write a pipeline for ChipSeq using Bowtie2. I'm trying to convert fastq file into sorted bam file, without duplicates. I wrote the following command:
${bowtie2_source} -x ${indexed_hg18} -U ${fastq_file} -S | ${samtools} view -bS | ${samtools} sort | ${samtools} rmdup -S ${target_dir}/${sample_name}.bam
Is it Ok in your opinion? or do I need to add at the end also -o ?
samtools view is not necessary, you can pipe directly to samtools sort, if your samtools version is reasonably recent.
Also, do you need to use hg18? It's quite old.
Which samtools version are you using?
Should be okay. I'm not sure if this would matter, but I tend to add in a - to tell samtools we're piping and has to read from stdin.
Oh, and maybe you should add -O BAM to samtools sort. Usually samtools sort is the end of my piping chain, but you continue with rmdup.
I think that up to samtools sort it works fine. The problem is with samtools rmdup, just not sure what... I ran the following command:
${bowtie2_source} -x ${ref_genome} -U ${fastq_file} -S | ${samtools} sort | ${samtools} rmdup -S - ${target_dir}/${sample_name}.bam
and received the error: [bam_rmdup] input SAM does not have header. Abort!
Then you probably need samtools sort -h
This question has been posted in a new thread so any discussion should be moved there for this specific question: Error with samtools rmdup
Use GATK: How to generate an unmapped BAM from FASTQ
Use BBMap: Reformat User Guide
EDIT: or use @Pierre's solution.
Log in to answer this question.
If an answer was helpful you should upvote it, if the answer resolved your question you should mark it as accepted. Please do the same for your previous posts as well.
Pierre has the right answer. The way you have written it, the software thinks that you want to look at all entries in the bam that have the chromosome name ${target_dir}/${sample_name}. This can't be done unless the sam/bam has an index, so the software is complaining.