This is a test version of Biostars. For the public version, visit https://www.biostars.org.
Bwa Error While Trimming Color Space Data

I am using the following pipeline for aligning some paired end SOLiD data. 

bwa aln -c -n 0.06 -o 2 -t 8 -q 10 ~/genomes/hydra/ACZUJGI/color/hydra ~/hydra/solid/hsamp_F3.fastq.gz > /scratch/hydra/hsamp_F3.sai
bwa aln -c -n 0.06 -o 2 -t 8 -q 10 ~/genomes/hydra/ACZUJGI/color/hydra ~/hydra/solid/hsamp_R3.fastq.gz > /scratch/hydra/hsamp_R3.sai

bwa sampe -P ~/genomes/hydra/ACZUJGI/color/hydra /scratch/hydra/hsamp_F3.sai /scratch/hydra/hsamp_R3.sai ~/hydra/solid/hsamp_F3.fastq.gz ~/hydra/solid/hsamp_R3.fastq.gz | samtools view -bS -|samtools sort - /scratch/hydra/hsamp_solid

On running this I get the following error. Just showing the last few lines from the output here. I only get this error when the -q parameter is nonzero. 

[bwa_paired_sw] 91 out of 33101 Q17 discordant pairs are fixed.
[bwa_sai2sam_pe_core] time elapses: 74.36 sec
[bwa_sai2sam_pe_core] refine gapped alignments... 1.53 sec
[bwa_sai2sam_pe_core] print alignments... [samopen] SAM header is present: 20914 sequences.
Parse error at line 20916: sequence and quality are inconsistent

The error happens at the conversion to bam step in the pipeline. If I look at line 20916, it shows

1_29_54 141     *       0       0       *       *       0       0       NNCANGNAANANATCNNCCGGNTANANTTGANTTANNTTN        !!@;!9!:?!;!:>>!!8?57!66!7!8=<9!<<9!!-?!!!&!<!!!!!      XC:i:40

Thus the read is truncated but the quality line is not. Is there a workaround for this?

bwa alignment
• 4,031 views • link
updated by brentp brentp • written by Farhat Farhat

1 answer

I think the proper answer is don't use -q with colorspace as it's designed for base-space. If you disregard that, you can pipe the output to this command (then to SAM)

$BWA_COMMAND \
| awk 'BEGIN{FS=OFS="\t"} \
      ($1 ~ /^@/){ print $0} \
      ($1 !~ /^@/){ $11 = substr($11, 0, length($10)); print $0}' \
| $SAMTOOLS_COMMAND > $OUT

That makes sure the qualities are the same length as the sequence.

You could also trim your reads with this: https://github.com/brentp/bio-playground/blob/master/solidstuff/solid-trimmer.py

• 1 views • link
updated by Ram Ram • written by brentp brentp

Thanks! I thought there may be some switch I am missing or something like that. The solid_trimmer looks useful!

• 0 views • link
written by Farhat Farhat

Log in to answer this question.

More posts like this
  • bwa mem vs bwa aln
    written by Adrian Pelin Adrian Pelin
    <p>Hello,</p> <pre> My data is 100bpx2 HiSeq PE reads, some PE overlap so I merged them into SE. For my analysis I used the aln …
  • Bwa Run With Fastq File With Different Read Length
    written by ttom ttom
    <p>I have 2 fastq files for the same sample, Sample1, but each with different read length(75bp and 100bp)</p> <p><strong><code>zcat Sample1_TTAGGC_L001_R1_001.fastq.gz | head -8</code></strong> # few …
  • Text Label In Circos
    written by dustar1986 dustar1986
    <p>Hi,</p> <p>I'm using Circos 0.62 to build up ideomgram. But I found that text label never showed (and "-debuggroup textplace" always returns "not_place" ).</p> <p>My …
  • Bwa-Mem And Sam-Flag "Read Paired In Proper Pair"
    written by Pierre Lindenbaum Pierre Lindenbaum
    <p>I'm playing with <strong>BWA-MEM</strong>.</p> <p>I've extracted a pair of reads properly paired (FLAG=2) with <a href='http://bio-bwa.sourceforge.net/'>bwa</a>-mem:</p> <pre><code>bwa mem human_g1k_v37.fasta NA12878_1.fastq.gz NA12878_2.fastq.gz | samtools view -S …
  • Bwa: Should Sai/Sam Sequences Processed Match?
    written by jomaco jomaco
    <p>Hi,</p> <p>After running <a href='http://bio-bwa.sourceforge.net/'>BWA</a> the output includes information such as this:</p> <pre><code>[bwa_sai2sam_pe_core] print alignments... 0.66 sec [bwa_sai2sam_pe_core] 197538574 sequences have been processed. [main] Version: …
  • Bwa Sampe Segmentation Fault
    written by a11msp a11msp
    <p>Hi everyone,</p> <p>I'm running <a href='http://bio-bwa.sourceforge.net/'>bwa</a> in the sampe mode and, after successfully processing &gt;10M reads it fails with a segmentation fault (as follows) on …
  • Using Bwa Sampe Using Illumna Pair End Reads
    written by Juliofdiaz Juliofdiaz
    <p>I am trying to map illumina paired reads using BWA and I am having some trouble understanding the workflow. I have two fastq files from …
  • Bwa - Bwasw Option - What'S The Output?
    written by Applet • Applet
    <p>I'm trying to give the user of my program the option to run <strong>bwa aln</strong> or <strong>bwa bwasw</strong>.</p> <p>But from reading the <a href='http://bio-bwa.sourceforge.net/bwa.shtml'>documentation</a> I …
  • Filtering Out Reads Based On Base Quality With Samtools Seems Not To Work Properly
    written by Pascal Pascal
    <p>We have an issue with the -Q parameter in samtools mpileup: it seems not working properly.</p> <p>In a bam file I know that the all …
  • Aligning Colorspace Reads Using Bwa
    written by Farhat Farhat
    <p>I have about 100 million reads from a SOLiD run. I am trying to align them using bwa and I got 0 alignments. What am …