Seqanswers Leaderboard Ad

Collapse

Announcement

Collapse
No announcement yet.
X
 
  • Filter
  • Time
  • Show
Clear All
new posts

  • Sequence and quality are inconsistent

    I hope this hasn't been addressed in a previous thread. I'm getting the following error message when I try to import into a .bam file (using samtools import).

    [samopen] SAM header is present: 84 sequences.
    Parse error at line 90: sequence and quality are inconsistent

    I'm using bwa aln to find coordinates and bwa sampe to generate alignments. When I look at that line I find the following record:

    HWI-ST609:107:C0DBJACXX:6:1101:2173:1948 83 2 112534540 0 66S35M = 112534381 -194 TTCACACCAGGGAGGAATTCCGAGTTCTTAAAGAGCCCCCGCTTAGGGTGGTTCTGCAGCCGCTCCTGATGGCTTCGGGAGCTGAAAAACTCTAACACNAA !
    ! RG:Z:UACC812BTR_tumor1 XC:i:35 XT:A:R NM:i:1 SM:i:0 AM:i:0 X0:i:2 X1:i:0 XM:i:1 XO:i:0 XG:i:0 MD:Z:32C2 XA:Z:2,-112534540,101M,1;


    I can go back to the fastq file and find the corresponding record:

    @HWI-ST609:107:C0DBJACXX:6:1101:2173:1948 1:Y:0:GCCAAT
    TTNGTGTTAGAGTTTTTCAGCTCCCGAAGCCATCAGGAGCGGCTGCAGAACCACCCTAAGCGGGGGCTCTTTAAGAACTCGGAATTCCTCCCTGGTGTGAA
    +
    ;;#2:5(@66:;;=<?8@)2<>>?@?39=??9>??3?>??;??2=?7>?999=9??#############################################

    Here are the aln, sampe and import commands I'm using:

    bwa-0.5.10/bwa aln -I -q 15 -t 12 human_g1k_v37.fasta XXXX_tumor1_1.fq.gz > XXXX_tumor1_read1.sai

    bwa-0.5.10/bwa aln -I -q 15 -t 12 human_g1k_v37.fasta XXXX_tumor1_2.fq.gz > XXXX_tumor1_read2.sai

    bwa-0.5.10/bwa sampe -P -r '@RG\tID:XXXX_tumor1\tLB:XXXX_tumor\tSM:XXXX_tumor' -f XXXX_tumor1.sam human_g1k_v37.fasta XXXX_tumor1_read1.sai XXXX_tumor1_read2.sai XXXX_tumor1_1.fq.gz XXXX_tumor1_2.fq.gz

    samtools-0.1.18/samtools import /human_g1k_v37.fasta.fai XXXX_tumor1.sam XXXX_tumor1_unsorted.bam

    I've tried using other versions of bwa and samtools. It seems the record is getting mixed up in the alignment process. Can anybody shed light on this?

    Myron Peto

  • #2
    The trivial answer...I think import is deprecated. So try SAMtools view; that command does a number of things, including converting sams to bams.

    Also, .sam files are really big. Where possible, you should pipe the output of whatever progam is making your .sam files (like bwa sampe, in your case) straight into samtools view.

    So your command line should look like

    bwa sampe ref.fa r1.sai r2.sai 1.fq 1.fq | samtools view -bS - > out.bam
    the -bS means "input is sam, output is bam".

    Comment


    • #3
      Good to know...

      Thanks for the response. I tried using samtools view on the data piped directly from bwq sampe and also on the sam file I generated. I still get the sequence and quality inconsistent error. From the record in the .sam file it seems the sequence and quality strings are of different length and the problem is upstream of the conversion from .sam to .bam

      Regardless, it's nice to know that view is the preferred option now.

      Myron

      Comment


      • #4
        Check to make sure index is compatible with 0.5.10 version of bwa. Indexing changed with 0.6.0.

        The previous aln or sampe command(s) is/are messing up.

        Your symptom is often the result of piping stderr to same same stdout in a previous process. Try adding " 2> stderr.txt" to explicitly force the stderr output to a file you can later examine on every command.

        Using a batch feeder that automagically adds a "2>&1" will pump stderr into stdout and you'll get that "mismatch read/quality thing " ... if your lucky.

        You can do a "strings filename" command on your sai file to see if there's anything bad there. ("bad"as in stderr messages, any other string should look like junk).

        Comment


        • #5
          Again, thanks for the response. I tried using 0.6.0 and 0.6.1 versions of bwa but got a segmentation fault with both of them. Maybe it would be useful to try to troubleshoot that path.

          As for the stderr and stdout, I left out the full version of the command. Here it is:

          bwa-0.5.10/bwa sampe -P -r '@RG\tID:XXXX_tumor1\tLB:XXXX_tumor\tSM:XXXX_tumor' -f XXXX_tumor1.sam human_g1k_v37.fasta XXXX_tumor1_read1.sai XXXX_tumor1_read2.sai XXXX_tumor1_1.fq.gz XXXX_tumor1_2.fq.gz 1>>XXXX_output.log 2>&1

          However, I've also tried running the aln and sampe commands without redirecting the outputs and I've gotten the same error every time.

          I tried strings on the .sai fill and did get a bunch of gibberish. I only went down a few pages and didn't grep the whole file for "error" or anything. I then piped it to "grep error" and came up empty.

          Myron

          Comment


          • #6
            Well. The -r parameter might be the problem. Bad output is RG:Z:UACC812BTR_tumor1 which looks suspiciously like '@RG\tID:XXXX_tumor1\tLB:XXXX_tumor\tSM:XXXX_tumor'

            Try running it without the -r and maybe -P parameters.
            It could be you don't have enough memory for -P.

            You can reheader the sam/bam output to put the RG in.

            Comment


            • #7
              Good suggestions they were, but, alas, I still get that same error message. I tried without the -r and them without either the -r or -P option but still get a misformed sam file.

              Any other ideas are most welcome.

              Myron

              Comment


              • #8
                were these files generated with CASAVA 1.8+?

                if so, try removing -I from command. This flag is for 1.5+ and now CASAVA can get qualities upto 41

                HTH

                dave


                #Edit: sorry, I meant 1.3+ and not 1.5+
                Last edited by dnusol; 02-03-2012, 04:57 AM.

                Comment


                • #9
                  Hallelujah!!

                  That last recommendation did the trick - the import (or view) from a sam to bam file is now proceeding. Knock on wood, hopefully it will finish without error.

                  Thanks much to all who offered suggestions, I offer my sincere gratitude.

                  Myron

                  Comment


                  • #10
                    For anyone who comes across this post in future with the same error message, one tip is to check the size of your .sai files. I had this error, one sai file was significantly smaller, I checked and one of the *.fastq.gz files had not quite transferred across correctly, hence the error.

                    Comment

                    Latest Articles

                    Collapse

                    • seqadmin
                      Exploring the Dynamics of the Tumor Microenvironment
                      by seqadmin




                      The complexity of cancer is clearly demonstrated in the diverse ecosystem of the tumor microenvironment (TME). The TME is made up of numerous cell types and its development begins with the changes that happen during oncogenesis. “Genomic mutations, copy number changes, epigenetic alterations, and alternative gene expression occur to varying degrees within the affected tumor cells,” explained Andrea O’Hara, Ph.D., Strategic Technical Specialist at Azenta. “As...
                      07-08-2024, 03:19 PM
                    • seqadmin
                      Exploring Human Diversity Through Large-Scale Omics
                      by seqadmin


                      In 2003, researchers from the Human Genome Project (HGP) announced the most comprehensive genome to date1. Although the genome wasn’t fully completed until nearly 20 years later2, numerous large-scale projects, such as the International HapMap Project and 1000 Genomes Project, continued the HGP's work, capturing extensive variation and genomic diversity within humans. Recently, newer initiatives have significantly increased in scale and expanded beyond genomics, offering a more detailed...
                      06-25-2024, 06:43 AM

                    ad_right_rmr

                    Collapse

                    News

                    Collapse

                    Topics Statistics Last Post
                    Started by seqadmin, 07-10-2024, 07:30 AM
                    0 responses
                    23 views
                    0 likes
                    Last Post seqadmin  
                    Started by seqadmin, 07-03-2024, 09:45 AM
                    0 responses
                    198 views
                    0 likes
                    Last Post seqadmin  
                    Started by seqadmin, 07-03-2024, 08:54 AM
                    0 responses
                    209 views
                    0 likes
                    Last Post seqadmin  
                    Started by seqadmin, 07-02-2024, 03:00 PM
                    0 responses
                    191 views
                    0 likes
                    Last Post seqadmin  
                    Working...
                    X